Avazuが終わった

Avazu Click-Through Rate Predictionが終わりました。長かった。期間が長すぎる!! と言っても、僕は始めて1週間目くらいにプランクトン分類が始まったため、そっちばかりでほとんどやっていなかったのですが、ついに終わりました。
結果は72th/1604。一応、上位5%には入っているけど、僕は絶対値を気にするので、72というのはあまり良くないな、と思う。

ソリューション

以下自分の投稿内容です。

概要

0.6 * 統計量を特徴量に使ったGBDTモデル + 0.4 * 交互作用項を追加した fast_solution_v3.py
です。

統計量を特徴量に使ったGBDTモデル

まずデータを2つに分割します。

  • 統計量計算用 (14/10/30のデータ)
  • モデル学習用 (14/10/30より前のデータ)

各カテゴリ変数について、

  • 出現回数
  • クリック回数 (出現したうちクリックされたレコード数)

を統計量計算用データから計算します。

カテゴリ変数 出現回数 クリック回数
site_id=123 10 3
site_id=456 100 10
app_id=123 40 8
app_id=456 50 10
... ... ...

という感じのデータ。これらを使って、モデル学習用データの各変数を置換します。
例えば、

site_id=123のとき site_id=(10, 3)
site_id=456のとき site_id=(100, 10)
app_id=123のとき app_id=(40, 8)
app_id=456のとき app_id=(50, 10)

のように。これでたくさんあるダミー変数の次元が全て2に次元削減されます。今回の場合、22変数だったので、44次元のデータになります。

この特徴量をxgboostを使って学習しました。
単体のLBスコアは、0.3933798でした。

今回とにかく、GBDTを使ってみたいという思いが初めにあったのですが、変数が多すぎてうまく学習できないように思えたのでとにかく減らしてみました。
データの分割方法はいろんなパターンが考えられて、たとえば、統計量計算に使うデータを1日ずつズラしながら7モデル学習して平均するようなこともできるのですが、やってみたところ特に良くなりませんでした。全体をランダムに分割するというのも試しましたが、これは統計量の情報が正確すぎるせいか、めちゃくちゃオーバーフィッティングしました。また古いデータから統計量を計算すると、新しく出て来た変数の情報が無くて変換できない変数がテストデータに増えるため、一番新しいデータである30日を統計量計算用に使うのがよさそうでした。

交互作用項を追加したfast_solution_v3.py

fast_solution_v3.pyはフォーラムで共有されているコードです。これ
中身は、FTRL-Proximal online learning algorithmとかいうのを使った線形分類器です。

変数の削除は、L1正則化に期待することにして、あとはいい感じのinteractionを追加できれば精度が良くなるだろうと思ったので、変数を追加する方向で改造しました。
方法としては、ただのブルートフォースアタックです。

  • 22変数から2つ選んだ全ての組み合わせについて、交互作用項として入れたときと入れなかった時の検証スコアを計算
  • 効果があった変数順にソート
  • 効果がありそうな変数を上から5個単位で追加していってスコアが悪くなったところで追加を辞める

というコードを書いて、一晩ほど計算した結果を少し手で調節して、以下のパラメーターで学習するとよい、ということになりました。

# これは事前に改造なしで探索したパラメーター
params = { "epoch": 2, "alpha": 0.05, "beta": 1.0, "L1": 3.0, "L2": 5.0 } 
# 追加する変数
params["interaction"] = (
    ('site_id', 'C14'),
    ('site_domain', 'C14'),
    ('app_domain', 'C14'),
    ('site_domain', 'C17'),
    ('site_id', 'C17'),
    ('app_domain', 'C17'),
    ('site_category', 'C14'),
    ('site_domain', 'C19'),
    ('app_id', 'C14'),
    ('site_domain', 'C20'),          
    ('app_id', 'C20'),
    ('site_id', 'C20'),
    ('site_id', 'C19'),
    ('app_category', 'C14'),
    ('app_domain', 'C18'),
    ('site_category', 'C17'),
    ('site_domain', 'app_id'),
    ('app_id', 'C17'),
    ('site_category', 'C20'),
    ('site_id', 'app_id'),
    ('site_domain', 'app_category'),
    ('app_id', 'C21'),
    ('site_domain', 'site_category'),
    ('app_id', 'C19'),
    ('site_category', 'app_id'),
    ('site_id', 'app_category'),
    ('site_domain', 'C21'),
    ('C14', 'C19'),
    ('site_id', 'site_domain')
    )

interactionは元のコードでは全て入れるか全て入れないかの2値ですが、指定した変数の組み合わせだけ使うように改造しています。
これで、単体のLBスコアは 0.3930048 でした。

アンサンブル

0.6 * 統計量を特徴量に使ったGBDTモデル + 0.4 * 交互作用項を追加したfast_solution_v3.py
です。重みは適当に試して決めました。

これで、LBスコアは 0.3905856 でした。各モデル単体のスコアに比べてかなり良くなっています。(このコンペではスコアが0.001変わると順位が150位くらい変わるので"かなり"良くなっています)

Privateスコアは 0.3887624。

感想

あまり熱心に取り組んでいないというのもあるのですが、フォーラムを見ると、learning rate decayが効いただの書いてあって、何でオレはそんなことも試していなんだ...と思う内容が多かったです。
データが大きくてあまり手の込んだことはできないという思いが強くて、試すパターンを制限しすぎてしまっていた気もします。

Torch7の分かりにくい話

Facebookがtorchの拡張を公開したこともあり、torchを使ってみている方も出始めたようなので、自分がtorchで使う上で分かりにくかったことやハマったことなどを上げていきます。
ある程度使っている方を前提としています。

ひどいことをたくさん書きますが、torchが嫌いなわけではありません。torchはサイコーです。

リリースやバージョンについて

ezinstallはgitのリポジトリをcloneしてインストールします。torchは今のところリリースの管理などしていなくて、日々更新される開発リポジトリがあるだけなので、インストールした日によって異なるバージョンがインストールされます。
たまに以下のようなコマンドを実行して更新したほうがいいです。

#!/bin/sh
sudo luarocks install cwrap
sudo luarocks install torch
sudo luarocks install optim
sudo luarocks install nn
sudo luarocks install image
sudo luarocks install cutorch
sudo luarocks install cunn
# 他にもありますが、自分はこのあたりだけ更新しています
# 更新順に気をつけないと依存関係が壊れることがあります

バージョンの管理が必要な場合は自分で行いましょう...

nnモジュール

学習用モデルと評価(予測)用モデルの切り替え

module:training() module:evaluate() で切り替えます。これは主にDropoutの動作に影響します。module:training()(デフォルトの状態)のままforward()で予測すると、Dropoutモジュールが含まれている場合にランダムに接続を切ってしまいます。予測を行う際は、必ず module:evaluate() でモードを切り替えてからforward()を行います。

バッチモードとはなにか

nnのほとんどのモジュール(全てではない..)はバッチモードを持っています。これはminibatch更新などで一度に複数のデータをforward/backwardするためのモードです。
このモードを使いたい場合は、入力として想定されているTensorよりひとつ軸の多いTensorを入力します。
たとえば、nn.Linear()は1D Tensorを想定しているモジュールなので1D Tensorを渡した時は逐次処理、2D Tensorを渡した場合はバッチ処理になります。nn.SpatialConvolutionMM()やnn.SpatialMaxPooling()は3D Tensorを想定しているので、3D Tensorを渡した時は逐次処理、4D Tensorを渡した時はバッチ処理になります。またCUDAバッグエンドのいくつかはバッチモードにしか対応していません。CUDAを使う場合はバッチモードで使うことを想定しているようです。

Criterion

Criterionは目的関数の定義です。が、SoftMaxなどの出力関数がCriterionと切り離されているため、分かりにくくなっていると思います。
SoftMax + MSEで最適化したい場合は、nn.SoftMax() + nn.MSECriterion() を使います。この場合、教師信号は、クラス数次元の1D Tensorになります(バッチモードの場合、2DTensor)。
SoftMax + NLLで最適化したい場合は、nn.LogSoftMax() + nn.ClassNLLCriterion()を使います。この場合、教師信号は、クラス番号(1からの整数値)になります(バッチモードの場合、1D Long Tensor)。また、nn.LogSoftMax()は、SoftMaxのlog、つまりexpする前の値を出力するため、確率を得たい場合は、forward():exp() で変換する必要があります。おそらく計算効率の問題でこうなっているのだと思いますがかなり分かりにくいと思います。
自分は、最近は、nn.SoftMax() + nn.TrueNLLCriterion()を使うようにしています。TrueNLLCriterionは、SoftMaxの出力をNegative Log Lossで直接最小化するというごく普通の動きをしてくれます。このモジュールは、facebookのfbnnの一部ですが、単体で動くのでこのファイルだけコピってきて使うのがオススメです。
https://github.com/facebook/fbnn/blob/master/fbnn/TrueNLLCriterion.lua

CUDA

CUDAを使う場合、Module, Criterionに対して:cuda()メソッドを呼び出して重み等のパラメーターをCudaTensorにするのと, 入力データを(:cuda()メソッドで)CudaTensorにして実行します。またforwardの返り値はCudaTensorになるのでその後使う場合はfloat()でFloatTensorに変換してから使ったほうがいいです。

例。

require 'cutorch'
require 'cunn'

local function build_model()
   -- 適当なCNN
   local model = nn.Sequential() 
   model:add(nn.SpatialConvolutionMM(3, 128, 5, 5, 1, 1))
   model:add(nn.ReLU())
   model:add(nn.SpatialMaxPooling(2, 2, 2, 2))
   model:add(nn.View(10 * 10 * 128))
   model:add(nn.Linear(10 * 10 * 128, 1024))
   model:add(nn.ReLU())
   model:add(nn.Dropout(0.5))
   model:add(nn.Linear(1024, 10))
   model:add(nn.LogSoftMax())
   
   return model
end
local function cuda_test()
   -- Module, Criterion, x, yをCudaTensorにする
   local model = build_model():cuda()
   local criterion = nn.ClassNLLCriterion():cuda()
   local x = torch.Tensor(64, 3, 24, 24):uniform():cuda() -- 入力
   local y = torch.Tensor(64):random(1, 10):cuda() -- ラベル
   -- batch forward
   local z = model:forward(x)
   local err = criterion:forward(z, y)
   -- batch backward
   model:backward(x, criterion:backward(z, y))
   print("CUDA Test Successful!")
end
torch.setdefaulttensortype('torch.FloatTensor')
cuda_test()
Dropout

Dropoutは第2引数にtrueを入れると学習時にランダムに接続を切って、予測時は(1-接続を切る確率)倍にするという一般的なDropoutになります。デフォルトでは、学習時にランダムに接続を切って1/(1-接続を切る確率)倍にして、予測時は出力そのままを使う実装になっています。デフォルトのほうが計算効率はいいと思うのですが、出力のスケールが変わるので他の実装を参考にした場合などは注意が必要です。

Convolution moduleいろいろ

SpatialConvolution/SpatialMaxPooling/SpatialAveragePoolingには沢山の種類があります。しかも微妙に動作が異なるものがあるので使う際には注意が必要です。

nn.SpatialConvolutionMM / nn.SpatialMaxPooling / nn.SpatialAveragePooling

torch標準のものです。他の実装に比べて遅いです(と言ってもすごく遅いわけではありません)。CPU,GPU,逐次,バッチと広く対応しています。
nn.SpatialAveragePoolingは平均ではなく合計を出力します。(平均を出力するものもあり、実装を切り替えただけつもりが出力が変わることがあります)
またこれらのモジュールは、出力サイズの計算にfloor(切り捨て)を使います。(ceil(切り上げ)の実装がいくつかあるため、実装切り替えただけつもりが出力のサイズが変わることがあります)

nn.SpatialConvolution

昔からあるものでサンプルなどで使われていることがありますが、最近はnn.SpatialConvolutionMMを使うようになっているので、使わないほうがいいです。

nn.SpatialConvolutionCUDA / nn.SpatialMaxPoolingCUDA

謎の実装です。おそらく実験的な実装のまま放置されているものなので使わないほうがいいです。

cudnn.SpatialConvolution / cudnn.SpatialMaxPooling / cudnn.SpatialAveragePooling

soumith/cudnn.torch · GitHub

NVIDIA cuDNNを使う実装です。luarocks install cudnnでインストールできます。別途NVIDIA cuDNNが必要です。NVIDIA cuDNNは最近公開されたRC2-v2以外は全てバグっているので必ず新しいバージョンを使うようにします。
実行は速いです。ただ、モデル保存時に100行くらいの警告が標準出力に出てきて、とてもうざいです。

ccn2.SpatialConvolution / ccn2.SpatialMaxPooling / ccn2.SpatialAvgPooling

soumith/cuda-convnet2.torch · GitHub

cuda-convnet2のカーネルを使う実装です。luarocks install ccn2でインストールできます。cuda-convnet2のコードは含まれているため別途インストールする必要はありません。おそらく一番速いです。(fbcunnは除く..)
これらのカーネルは、torch標準とはTensorのフォーマットが違います。torch標準では(batch, depth, height, width)ですが、ccn2は(depth, height, width, batch)です。nn.Transpose()で軸を入れ替えてから使います。
CUDA&バッチモードにしか対応していません。またフィルタの数は16の倍数、バッチの数は32の倍数、などの制限があります。
自分は主にこれを使っています。

nn.SpatialConvolutionCuFFT

facebookが公開した実装です。fbcunnに入っています。
カーネルが大きい場合は、とても速いらしいですが、小さなカーネル(流行りの3x3など)の場合はcudnnやcuda-convnet2のほうが速いです。またメモリ使用量が多いためTitan BlackやK40などメモリが多いデバイスでないと大きなモデルが学習できないことがあるようです。

他にもありますが、上のやつだけ把握していればいいように思います。
soumith/convnet-benchmarks · GitHub

に各実装のベンチマークがあります。

画像関係

image.loadで読み込んだ画像の各画素は0.0-1.0にスケールされています。0-255ではありません。

image.display()

image.display()で画像が表示できると書いてありますが、thコマンドではできません。
qluaコマンドを使うとqtが使えるようになっていて、image.display()などが使えます。
ただqluaはluajitと異なるバナーが表示され、luaの実装が違う?ようなので学習など重要な処理には使わないほうがいいと思います。
またimage.display()は明るさを勝手に調節するため、画像処理結果の確認はimage.save()で保存してから見たほうがいいです。

image.crop

torch.Tensorの添字は1オリジンですが、どうもimageモジュールの一部は0オリジンで処理されているところがあるようです。
image.cropもx/yは0から始まっているように見えます。

image.rotate

thetaの符号が逆?
領域外は黒(0)で塗られます。data augmentationなどで使う場合は、可能なら背景を0にするか、不可能なら広めに取った領域を回転させたあと領域外を含まない範囲でcropするのがいいと思います。
image.warpを使うと、領域外は端のピクセルを引き伸ばすようにできますが、変なパターンが入るのは同じことなので...

image.warp

ソース画像の変換後の座標を(2(x,y), height, width)のTensorで定義して画像を変換する関数です。大概の変換処理はこの関数でできますが、変換用のデータを自分で計算する必要があります。
いまのところbicubicに対応しているのはこの関数だけなので、rotateやscaleを使わずこの関数を使っている方が多いようです。
使い方は、https://github.com/torch/image/blob/master/test/test_warp.lua

これも0オリジンだと思います。(テストは1オリジンで計算していますが端のピクセルが領域外になっている)

gfx.js

一部のチュートリアルでは gfx.js を使って画像の表示等を行うと書いてありますが、現在メンテされておらずバグっていて動かないという話なので無視しましょう。開発側としてはiTorchに移行したいようです。

optim

config

optim.sgd()にlearningRateなどの設定を渡しますが、このテーブルはただ設定を受け取るだけはなく、次の更新時に使うパラメーター等が書き込まれます。なので

optim.sgd(a, b, {learningRate = 0.1, momentum = 0.9})

のようにリテラルや一時変数で渡してしまうと、momentumやlearningRateDecayの計算ができなくなってしまいます。学習中生存している変数を介して渡すようにします。

NAG

optim.sgdにnesterovというオプションがありこれをtrueにするとNesterov momentumになると書いてありますが、一般的な実装とは違うようです。この問題が上がったあとoptim.nagが追加されたのでNAGを使いたい場合はこちらを使ったほうがいいようです。

torch.setdefaulttensortype

一部のチュートリアルでは、CUDAを使う場合はデフォルトをCudaTensorにすると書いてありますが、それをやるといろいろなところがバグるためしてはいけません。
CudaTensorはまだFloatTensorと完全な互換性がありません。努力はしているようです。
missing API against TH · Issue #70 · torch/cutorch · GitHub
デフォルトはDoubleTensorですが、Floatで十分なので(どうせCUDAはFloatしか使わない)、FloatTensorにするのがオススメです。必要なメモリが半分になり処理も速くなります。

torch.manualSeed

CUDAを使っているときは、cutorch.manualSeedも設定します。

便利情報

torchサイアクだという印象を残して終わらたないために、最後に便利情報を載せます。

Caffe関連

Caffeうらやましい!と思うところはだいたい取り込まれています。

szagoruyko/loadcaffe · GitHub
を使うと、Caffeのpretrained imagenet modelをtorchのモジュールに変換して読み込めます。(caffeのインストールは必要ありません)

szagoruyko/torch-caffe-binding · GitHub
はcaffeをtorchのモジュールとして使います(caffeのインストールが必要です)。caffeはtorchの部品のひとつであると言えます。

画像認識関連

画像認識は、ImageNet Challengeの成果が非常に参考になります。

convnet-benchmarks/torch7/imagenet_winners at master · soumith/convnet-benchmarks · GitHub
には、AlexNet[1], Overfeat[2], GoogLeNet(補助分類器なし)[3], VGG model[4]の実装例があります。(速度のベンチマーク用なので完全でないものがあります)

ImageNet-Training/Models at master · eladhoffer/ImageNet-Training · GitHub
には、AlexNet, GoogLeNet, NIN[5], OverFeatの実装例があります。

ILSVRC2014の資料は
ImageNet Large Scale Visual Recognition Competition 2014 (ILSVRC2014) Workshop
にあります。

あとオレのCIFAR-10のコードも参考になります。(オチです)
Kaggle CIFAR-10の話 - デー

Torch7で最近傍探索のベンチマーク

最近、Torch7で最近傍探索を繰り返し行いたかったけど、すごく遅いのでは??という不安があったのでk-NN(k=1)でベンチマークしてみた。

設定

  • MNISTをk-NN(k=1)で評価する
  • 尺度はコサイン類似度とする
  • テストを全件評価してかかった時間を計測する

環境

  • Intel(R) Core(TM) i7-3770K CPU @ 3.50GHz
  • 32GB RAM
  • GeForce GTX 760

パッと思いついた実装

初めは難しいことは考えず、パッと思いついた方法を試してみる。

工夫としては、

  • コサイン類似度を求める際にベクトルのノルムを毎回計算したくないので、最初にノルムが1になるように正規化しておく(内積=コサイン類似度になる)
  • 各テストデータと各学習データの比較は、gemv(torch.mv)で一度に計算すれば速いのではないか
require 'optim'
require 'xlua'

torch.setdefaulttensortype('torch.FloatTensor')
local EPSILON = 1.0e-6

-- normを1に正規化
local function normalize_l2(x)
   local norm = torch.pow(x, 2):sum(2):sqrt():add(EPSILON)
   x:cdiv(torch.expand(norm, x:size(1), x:size(2)))
end

function main()
   local mnist = require 'mnist'
   local trainset = mnist.traindataset()
   local testset = mnist.testdataset()
   local train_x, train_y = trainset.data:float(), trainset.label
   local test_x, test_y = testset.data:float(), testset.label
   local classes = {0,1,2,3,4,5,6,7,8,9}
   local confusion = optim.ConfusionMatrix(classes)
   local t = sys.clock()
   
   -- データを整形
   train_x = train_x:reshape(train_x:size(1), 28 * 28)
   test_x = test_x:reshape(test_x:size(1), 28 * 28)
   train_y:add(1)
   test_y:add(1)
   
   -- L2 normを1に正規化
   normalize_l2(train_x)
   normalize_l2(test_x)
   
   -- 各テストデータについて
   for i = 1, test_x:size(1) do
      -- コサイン類似度が最も大きいインスタンスを選択
      local _, nn_index = torch.mv(train_x, test_x[i]):max(1)
      -- 結果を評価
      local y = train_y[nn_index[1]]
      confusion:add(y, test_y[i])
      if i % 100 == 0 then
	 xlua.progress(i, test_x:size(1))
      end
   end
   -- 結果を表示
   print(confusion)
   print(string.format("*** %.2fs", sys.clock() - t))
end
main()
結果
ConfusionMatrix:
[[     978       1       0       0       0       0       0       1       0       0]   99.796% 	[class: 0]
 [       0    1129       3       1       0       1       1       0       0       0]   99.471% 	[class: 1]
 [       9       0    1003       4       0       0       2      10       3       1]   97.190% 	[class: 2]
 [       0       0       1     977       0      13       0       5       9       5]   96.733% 	[class: 3]
 [       1       3       0       0     940       0       6       3       1      28]   95.723% 	[class: 4]
 [       1       1       0      17       1     852      10       1       4       5]   95.516% 	[class: 5]
 [       4       3       0       0       2       3     946       0       0       0]   98.747% 	[class: 6]
 [       2      11       5       2       2       0       0     995       0      11]   96.790% 	[class: 7]
 [       6       1       1      13       2       3       5       4     935       4]   95.996% 	[class: 8]
 [       5       6       1       4       9       3       1       8       4     968]]  95.937% 	[class: 9]
 + average row correct: 97.189832925797% 
 + average rowUcol correct (VOC measure): 94.583150148392% 
 + global correct: 97.23%
*** 95.11s	

95.11秒だった。正解率は97.23%。世の中にはMNISTをk-NNしてみたら数時間かかったとか言っている人も散見されるので、そう考えると速い気もする。

gemvではなくgemmで計算する

インスタンスひとつづずgemvしてたけど、全体をgemm(torch.mm, 行列の積)で計算すればもっと速くなるだろうと思ったので変更してみた。

require 'optim'
require 'xlua'

torch.setdefaulttensortype('torch.FloatTensor')
local EPSILON = 1.0e-6

-- normを1に正規化
local function normalize_l2(x)
   local norm = torch.pow(x, 2):sum(2):sqrt()
   norm:add(EPSILON)
   x:cdiv(torch.expand(norm, x:size(1), x:size(2)))
end

function main()
   local mnist = require 'mnist'
   local trainset = mnist.traindataset()
   local testset = mnist.testdataset()
   local train_x, train_y = trainset.data:float(), trainset.label
   local test_x, test_y = testset.data:float(), testset.label
   local classes = {0,1,2,3,4,5,6,7,8,9}
   local confusion = optim.ConfusionMatrix(classes)
   local t = sys.clock()
   
   -- データを整形
   train_x = train_x:reshape(train_x:size(1), 28 * 28)
   test_x = test_x:reshape(test_x:size(1), 28 * 28)
   train_y:add(1)
   test_y:add(1)
   
   -- L2 normを1に正規化
   normalize_l2(train_x)
   normalize_l2(test_x)
   -- 全てのテストデータについてコサイン類似度が最も大きいインスタンスを選択
   local cosine = torch.mm(train_x, test_x:t())
   local _, nn_index = cosine:max(1)
   -- 結果を評価
   for i = 1, test_x:size(1) do
      local y = train_y[nn_index[1][i]]
      confusion:add(y, test_y[i])
      if i % 100 == 0 then
	 xlua.progress(i, test_x:size(1))
      end
   end
   -- 結果を表示
   print(confusion)
   print(string.format("*** %.2fs", sys.clock() - t))
end
main()
結果
ConfusionMatrix:
[[     978       1       0       0       0       0       0       1       0       0]   99.796% 	[class: 0]
 [       0    1129       3       1       0       1       1       0       0       0]   99.471% 	[class: 1]
 [       9       0    1003       4       0       0       2      10       3       1]   97.190% 	[class: 2]
 [       0       0       1     977       0      13       0       5       9       5]   96.733% 	[class: 3]
 [       1       3       0       0     940       0       6       3       1      28]   95.723% 	[class: 4]
 [       1       1       0      17       1     852      10       1       4       5]   95.516% 	[class: 5]
 [       4       3       0       0       2       3     946       0       0       0]   98.747% 	[class: 6]
 [       2      11       5       2       2       0       0     995       0      11]   96.790% 	[class: 7]
 [       6       1       1      13       2       3       5       4     935       4]   95.996% 	[class: 8]
 [       5       6       1       4       9       3       1       8       4     968]]  95.937% 	[class: 9]
 + average row correct: 97.189832925797% 
 + average rowUcol correct (VOC measure): 94.583150148392% 
 + global correct: 97.23%
*** 12.41s

12.41秒だった。正解率は当然同じ。かなり速くなった。

CUDAでやってみる

Torch7はBLASでできるような計算なら簡単にCUDA化できるのでやってみた。

工夫として、

  • gemmは使用メモリが多すぎてGPUにメモリが確保できなかったので分割して計算するようにした
require 'cutorch'
require 'optim'
require 'xlua'

torch.setdefaulttensortype('torch.FloatTensor')
local EPSILON = 1.0e-6

-- normを1に正規化
local function normalize_l2(x)
   local norm = torch.pow(x, 2):sum(2):sqrt()
   norm:add(EPSILON)
   x:cdiv(torch.expand(norm, x:size(1), x:size(2)))
end

-- 行列の積を直接計算しようとするとGPUメモリに載らなかったので16分割して計算
local function split_mm(a, b)
   local BLOCKS = 16 -- 分割数
   local step = math.floor(b:size(1) / BLOCKS)
   local results = torch.Tensor(a:size(1), b:size(1))
   for i = 1, b:size(1), step do
      local n = step
      if i + n > b:size(1) then
	 n = b:size(1) - i
      end
      if n > 0 then
	 results:narrow(2, i, n):copy(torch.mm(a, b:narrow(1, i, n):t()))
      end
      collectgarbage()
   end
   return results
end

function main()
   local mnist = require 'mnist'
   local trainset = mnist.traindataset()
   local testset = mnist.testdataset()
   local train_x, train_y = trainset.data:float(), trainset.label
   local test_x, test_y = testset.data:float(), testset.label
   local classes = {0,1,2,3,4,5,6,7,8,9}
   local confusion = optim.ConfusionMatrix(classes)
   local t = sys.clock()

   -- データを整形
   train_x = train_x:reshape(train_x:size(1), 28 * 28)
   test_x = test_x:reshape(test_x:size(1), 28 * 28)
   train_y:add(1)
   test_y:add(1)
   
   -- L2 normを1に正規化
   normalize_l2(train_x)
   normalize_l2(test_x)
   
   -- 計算用のデータをCudaTensorに変換(GPUのデバイスメモリに転送)
   train_x = train_x:cuda()
   test_x = test_x:cuda()
   
   -- 全てのテストデータについてコサイン類似度が最も大きいインスタンスを選択
   local cosine = split_mm(train_x, test_x)
   local _, nn_index = cosine:max(1)
   -- 結果を評価
   for i = 1, test_x:size(1) do
      local y = train_y[nn_index[1][i]]
      confusion:add(y, test_y[i])
      if i % 100 == 0 then
	 xlua.progress(i, test_x:size(1))
      end
   end
   print(confusion)
   print(string.format("*** %.2fs", sys.clock() - t))
end
main()
結果
ConfusionMatrix:
[[     978       1       0       0       0       0       0       1       0       0]   99.796% 	[class: 0]
 [       0    1129       3       1       0       1       1       0       0       0]   99.471% 	[class: 1]
 [       9       0    1003       4       0       0       2      10       3       1]   97.190% 	[class: 2]
 [       0       0       1     977       0      13       0       5       9       5]   96.733% 	[class: 3]
 [       1       3       0       0     940       0       6       3       1      28]   95.723% 	[class: 4]
 [       1       1       0      17       1     852      10       1       4       5]   95.516% 	[class: 5]
 [       4       3       0       0       2       4     945       0       0       0]   98.643% 	[class: 6]
 [       2      11       5       2       2       0       0     995       0      11]   96.790% 	[class: 7]
 [       6       1       1      13       2       3       5       4     935       4]   95.996% 	[class: 8]
 [       5       6       1       4       9       3       1       8       4     968]]  95.937% 	[class: 9]
 + average row correct: 97.179394364357% 
 + average rowUcol correct (VOC measure): 94.562811851501% 
 + global correct: 97.22%
*** 5.12s	

5.12秒だった。爆速!

プランクトンの画像分類コンペが開催されています!

Kaggleでプランクトンの画像分類モデルの良さを競うコンペが開催されています。
Description - National Data Science Bowl | Kaggle

Computer VisionやDeep Learningに興味のある方、自分の力を試してみたい方、手頃なサンドバックが欲しい方などは参加されてみてはいかがでしょうか。

概要(超訳)
海の状態をプランクトンの種類や数の分布などから調べたい。プランクトンの識別を人間が行うのは大変なので水中カメラで撮影したプランクトンの画像から種類を予測するモデルを作ってもらいたい。
期間
2014/12/15〜2015/3/16
賞金総額
$175,000 (14/12/20: 約2090万円)
データ
学習画像: 約3万枚 (グレースケール/様々なサイズ/不均一), 対象クラス: 121種類(多クラス分類/確率出力/multi-class logarithmic lossで評価), テスト画像: 約13万枚 (チート対策で評価に使用されないデータもたくさん含む)

詳細は、サイトを参照ください。
注意することとして、テストデータを使っての半教師あり学習や教師なし特徴学習は許可されていますが、外部データを使うことは禁止されています。(Caffeのpre-trained ImageNet modelなどは使用不可)

賞金が多いことから、たくさん人が来ると思いきや、まだあまりいないようなので、宣伝を書きました。当然僕も参加しています。

Kaggle CIFAR-10の話

以前、Kaggle CIFAR-10 に参加していると書きましたが、これが2週間ほど前に終わりました。コンペはまだ Validating Final Results の状態なのですが、2週間たっても終わらず、いつ終わるのか謎なのと、多分結果は変わらないと思うので先に書きます。

CIFAR-10は、次のような32x32の小さな画像にネコ、犬、鳥など10種類の物体が写っているので、与えられた画像に何が写っているか当てる問題です。
f:id:ultraist:20141108185409p:plain
(Kaggle CIFAR-10のデータセットは、通常のCIFAR-10と結果の互換性がありますが、チート防止に画像のハッシュ値が変わるように改変されているのと、テストセットに29万枚のジャンクイメージが含まれています。)

自分の結果は、0.9415 (正解率94.15%)で、Classification datasets results によると、state-of-the-artが91.78%なので、それを上回って人間による精度である94%に達しているのですが、このスコアでなんと5位でした。1位はDeepCNetで、95.53%という驚異の精度を出しています。2位もDeepCNetとDropConnectの結果を合わせたものなので、DeepCNet最強だったという感じです(DeepCNetのコードはコンペ終了前に公開されていました)。3、4位は手法を公開していないので不明です。

自分の手法

kaggle-cifar10-torch7 - Github
でコードを公開しています。Torch7で実装しています。

オリジナル性が高いものはなく、よくある手法をいくつか組み合わせたのと、VGG(University of OxfordのVisual Geometry Group)がILSVRC2014で使ったモデルをCIFAR-10に調節しただけものです。

やったことは

  • 学習データを36倍に増化(Data Augmentation)
  • GCN + ZCA Whiteningで正規化
  • VGGのモデルをベースにしたConvolutional Neural Network(CNN)を学習
  • 上記のモデルを重みの初期値とMini-Batch-SGDの更新順を変えて6個学習し、各分類器の平均を予測として出力

です。

学習データを36倍に増化(Data Augmentation)

ニューラルネットワークは、経験的にはオーバーフィッティングしなければ層数や素子数が多いほどよいというのがあって、よりDeeeepしたいという思いはありますがデータが少ないとオーバーフィッティングしてしまうので、データを増して複雑なモデルでもオーバーフィッティングしにくいようにしました。
コツとしては、できるだけ"あり得る範囲"の変換により増やすということです。画像の場合、人工的なノイズや歪めたりでいくらでもパターンが増やせますが、テストセットに出てこないようなパターンで学習データの分布を歪めてしまうとよくないので、元のデータとは違うけどテストには出てくるパターンに変換できるのが理想です。
今回は次の3つのメソッドを使いました。

Cropping
CIFAR-10の学習画像は32x32ですが、これを24x24の部分画像に分解します。4px置きに切り出すと、3x3の9パターンが切り出せるのでデータを9倍に増やせます。(学習画像は小さくなります)
Scaling
Croppingでの切り出しサイズを28x28にして2px置きに切り出すと3x3の9パターンが切り出せます。この部分画像を24x24にリサイズ(ズームアウト)して学習画像とします。
Horizontal reflection
左右反転です。これまで増やした画像を左右逆の2のパターンに分けて2倍に増やします。

これで(9 + 9) * 2 = 36倍になります。学習画像は32x32ではなく、24x24になります。
予測時は、予測対象の画像を同じ方法で36倍に増やして、各画像に対して予測を行い、それらを平均して予測結果としています。当然、予測にかかる時間も増えます。
この処理によって、2〜3%くらい精度がよくなりました。

f:id:ultraist:20141108184751p:plain
1行目がCropping、2行目がCropping+Scaling、3、4行目がHorizontal reflectionです。

GCN + ZCA Whiteningで正規化

GCN(Global Contrast Normalization)は、standardizeとかz-scoreとか言われるものと同じで、データ全体から各要素の平均と標準偏差を求めて、平均を引いて標準偏差で割ります。入力の値域が-2〜2くらいに正規化されて、スケールの異なる軸があった場合でもその範囲にそろえられます。また、よく出る値は平均に近くなり、あまり出ない値は大きな絶対値を持つようになります。スケールの大きな軸の影響を抑えるのと、学習時の収束が速くなる効果があります。
2014/12/20追記
この説明は間違っていました。GCNは、画像をまたがずに、画像内での平均と分散を求めて平均を引いて分散で割るようです。これは、"An Analysis of Single-Layer Networks"の実装では、local contrast normalizationと書かれていた処理で、z-scoreはglobal standardizationと書かれていて、Maxoutの論文では、このlocal contrast normalizationがglobal contrast normalizationと書かれていたので、globalとlocalの概念がどこにあるのか混乱して勘違いしていました。ただ自分の実装では、このブログ通りのz-scoreを使っています。

ZCA Whiteningは、データ全体の分散共分散行列の固有ベクトルで主軸変換を行なって、変換後の空間でstandardizeを行なって元の空間に戻すというものです。自然画像は、あるピクセルはその近隣のピクセルと相関が強いという特徴があるので、この相関を消すことで色の情報を持ったままエッジ検出を行ったような結果が得られます。元の空間に戻すのは、CNNが元画像の構造を前提としているからです。
ZCA Whiteningは、An Analysis of Single-Layer Networks in Unsupervised Feature Learning - Andrew Ngですごくよい結果を出した前処理で、僕もこの手法を実装したことがあるのですが、この手法においてはZCA Whiteningをするかしないかで、CIFAR-10の精度が15%くらい変わります。ただ、最近のDeep CNNではほとんど差がでないので、必要なかったのではないかと思っています。その前のモデル(Network In Network)ではほんの少しだけ精度が改善できていたのと、外して変わらない精度が出るか試している余裕がなかったので、そのまま慣性で入れています。

VGGのモデルをベースにしたDeep Convolutional Neural Networkを学習

[1409.1556] Very Deep Convolutional Networks for Large-Scale Image Recognition で提案されているものをベースにしたDeep Convolutional Neural Networkで分類器を作りました。

伝統的なCNNでCIFAR-10用のアーキテクチャを作ると、conv 5x5 -> maxpool -> conv 5x5 -> maxpool -> conv 5x5 -> maxpool -> fc(fully connected) -> softmaxのようになるのですが、このconv 5x5の部分を3x3カーネルを2つか3つ並べたものに置き換えます。
これで

  • 層数が増える
  • 線形の大きなカーネル非線形(convごとにReLUを挟んでいるため)の3x3を並べたものに置き換えるので表現力が上がる
  • 大きなカーネルで一回畳み込むよりも小さなカーネル複数回畳み込んだほうが計算量が少ない(5x5 > 3x3x2, 7x7 > 3x3x3)

というような効果があります。また3x3カーネルに1pxのpaddingを加えると、畳み込み層によって画像サイズが縮小されなくなるので、理論上は無限に層数を増やせるようになります。これはCIFAR-10のような入力画像が小さい場合に嬉しいです(畳込みでサイズが減っていくと増やせる層数に限界があるため)。

最終的に使ったアーキテクチャは、

conv 3x3 -> conv 3x3 -> maxpool -> conv 3x3 -> conv 3x3 -> maxpool -> conv 3x3 -> conv 3x3 -> conv 3x3 -> conv 3x3 -> maxpool -> fc -> fc -> softmax

というDeeeepなものです。詳しくはソースコードのページに表を書いているので参照してください。

上記のモデルを重みの初期値とMini-Batch-SGDの更新順を変えて6個学習し、各識別器の平均を予測として出力

ニューラルネットワークを使ったモデルで簡単に精度を上げる方法として、いくつかのモデルを学習して平均を取るというのがあります。ニューラルネットワークは初期値依存があるのと大域最適化はされないので、乱数のseedが異なると(微妙に)異なる結果を出力するモデルが学習されます。なので、いくつか学習して平均を取ると結果が安定します。

よくやるのはBagging(Committee Network)ですが、Baggingはサンプリングの割合など調節しなければならないのと今回はこれをやろうと思ったのが終了3日前で、調節する余裕がなく一発勝負だったため、良くなることはあっても悪くなることはないだろうという考えで以下の設定で行いました。

  • 同じ学習データ
  • 異なる初期重み
  • 異なる更新順

また今回使ったモデルは学習に20時間程度かかり、単体マシンで学習していては2つしか学習できないので、EC2のSpot InstanceでGPU Instanceを6個たち上げて並行して学習しました。

結果的には、シングルモデルだと93.33%、6モデルの平均で94.15%だったので、この処理で0.85%改善できていました。

その他の話

VGGの論文が発表される前は、Network In Networkを使っていました。これは最終的には92.4%の精度を出せたので、悪くはなかったと思います。

最後の方では、このままではどうやってもDeepCNetに勝てないと思ったので、GoogLeNetを実装してみたのですが、学習がクッソ遅い上validationで88%しか出なかったので諦めました(調節が足りないのか、問題に合っていないのか、なにか間違っているのか分かっていない)。

この2つのモデルは、参考としてソースコードのディレクトリに置いてあります。(nin_model.luaとinception_model.lua)

lbpcascade_animefaceをgithubに置きました

見る前からスター !lbpcascade_animeface · GitHub

OpenCV用のアニメ顔検出器 lbpcascade_animeface をgithubに置きました。
内容は2011年から変わっていません。
妙なブログ記事ではなく参照しやすくなったのではないかと思います。
全宇宙に拡散するためスターをお願いします。
ちなにみこれはImage::AnimeFaceとは全く別のものです。

見た後にもスター !lbpcascade_animeface · GitHub

AfSIS 反省会

Kaggleで行われていたAfrica Soil Property Prediction Channgleに参加して、結果は1233チーム中148位だった。
反省会です!!

ちなみに参加していなかった人には何を言っているか全く分からないと思いますのでご了承ください。

今回のコンペは、yag_aysさんも書いているけど、public LBとprivate LBで、かなり大きな順位変動があって、public LB1位だった人が500位台まで飛んだり、上位を狙うのは極めて難しい内容だった。ただこれは結果が出る前から明らかだったことで、それを見越してprivate LBで上位に入るように考えていたんだけど、いろいろ予想が外れてしまい失敗してしまった。

オレの予想と現実の違い

予想
  • private 1位は0.46くらい
  • (この間に多くても20人くらい)
  • オレは0.48〜0.50くらい
  • Abhishekのコードは0.55以上
現実
  • private 1位は0.46892
  • (この間に150人くらい)
  • オレは0.50290
  • Abhishekのコードは0.50558

現実は厳しい!!

だいたい以下のようなことを考えていた。

  • 自分の投稿やフォーラムの話、Sentinel LandscapeベースのCVからpublic LBより0.1くらい悪いスコアになるだろう
  • public LB上位(0.36~0.38)のうち何人かは本当にうまくいっていて残るだろうから1位は0.46くらい
  • 自分の投稿は、Sentinel LandscapeベースのCVで0.48くらい、public LBの計算に使われている簡単なデータセットはprivateから除かれるのでこれよりちょっと悪くなるだろう
  • AbhishekのコードはSentinel LandscapeベースのCVで0.61、locationベースのCVで0.56くらいという話があったので(自分では検証していない!!)よくても0.55くらいだろう
  • 参加者の多くは、Abhishekのコードをベースにしているので0.5には届かないだろう

多くの人はprivate LBで爆死するので、爆死しないように気をつけておけばそれほど頑張らなくても上位に入れるだろうと考えていたんだけど、この予想が外れ、Abhishekのコードがなんと0.50というスコアを持っていたため、そこから少しでも改善できた人達150人くらいに抜かれてしまった。

反省点としては、他の人が勝手に爆死することを期待したり、Kaggleの天才データサイエンティスト達を舐めてたりといったふざけた態度でのぞまず、自分のスコアを上げる方向で頑張るべきだった。マル。

自分のソリューション

ちなみに自分の方法は以下のような感じだった。

  1. spectra特徴からCO2の領域を削除
  2. spectraにfirst derivative filterを適用
  3. spectraをPCA Whiteningで160次元に圧縮
  4. 圧縮したspectra特徴, spatial特徴(衛星から取ったやつ), depth(表層か下層か)を入力可能な特徴とする
  5. huber loss functionのneural networkで回帰
  6. neural networkのアーキテクチャやハイパーパラメーターはSentinel Landscapeに基づくCVで良い結果が出るように調節する (spatial特徴を使うか使わないか、depthを使うか使わないかの組み合わせも含めて)

1は、ホストが推薦していたことで、やったほうが結果もよかったのでやった。
2は、CVではやってもやらなくてもほぼ違いはなかったので迷ったけど、データをプロットしてみるとデータごとに異なる謎の直流成分(定数のズレ)が入っていて、どういう理由で入ったのか分からないけど、おそらくホスト側でデータを正規化したときに入ったんじゃないかと思ったのと、ホストが用意しているベースラインのベンチマークコードではfirst derivativeを使っていて、これを使うと直流成分は消えるので、消したほうが安心だろうと思って適用することにした。
3もやるか迷ったものだけど、正規化のパターンとして

  • z-socre
  • ZCA Whitening
  • PCA Whitening (512,256,160,128 component選択)

を試してみて、z-socreは一番結果が悪かったので無しで、ZCAとPCAはほんのちょっとだけZCAのほうがよかった。ただ、今回のデータは学習データが約1000件で入力が約8000次元あるので、できるなら圧縮したほうがいいだろうと思っていたのと、入力次元を少なくするとneural networkの学習時間がかなり短くなので、実験回数が増やせて、この微かな差は挽回できるのではないか、と考えて160次元に圧縮することにした。
結果からするとこれは失敗だった。ZCAバージョンも投稿していたんだけど、private LBではZCAバージョンのほうが微かによくて(0.4991)、こっちを選んでいれば95位くらいだった。(そんな変わらない)

4.はspectra特徴だけで学習するか、spatialやdepthも含めるかという話だけど、これは目的変数ごとに検証した結果以下のようにした。
Ca: spectra
P: spectra + depth
pH: spectra + spatial + depth
SoC: spectra + spatial + depth
Sand: spectra + depth

5.は、この問題には非線形回帰がよいことはちょっとやってみれば分かることで、非線形回帰ができるモデルとして自分の得意なneural networkを選んだ。重要なのは損失関数で、MSEではなくHuber lossを使うことにした。これは、学習データの目的変数に明らかに大きすぎる外れ値っぽいデータが少しだけ入っていて、MSEを最小化するとその影響がすごく出てしまうので、ロバストにしたほうがいいだろうということで選んだ。CVでもMSEよりHuber lossのほうが良い結果が出ていた。
ちなみにSupport Vector Regression(Abhishekのコード)が良い結果を出していた理由も、損失関数にsquared lossを使っていないことだと思っていて、というか、SVRで良い結果が出るという話が出たあとで、なぜSVRだけ良い結果が出るのか考えた結果、ε-sensitive lossのおかげだろうと思ったので、自分もそれに習って性質が似ていてneural networkで扱いやすそうなHuber lossを使ってみることを思いついた。

6.は一番重要なところで、CVのfoldをランダムサンプリングではなく、Sentinel Landscapeという単位を壊さないように分けるようにする。
これが重要な理由は、データページに書いてある。

The training and test data have been split along Sentinel Landscape levels because we are primarily interested in predicting soil properties at new Sentinel Landscapes.

http://www.kaggle.com/c/afsis-soil-properties/data

学習データとテストデータは、Sentinel Landscape levelsで分けられていると明確に書いてある。Sentinel Landscapeは土壌をサンプリングした時の一番大きい空間単位で、

60 Sentinel Landscapes
16 Sampling Clusters per Sentinel Landscape
10 Sampling Plots per Sampling Cluster
2 composite Soil Samples (topsoil & subsoil) per Sampling Plot

と書かれている。つまりこのコンペのデータは、階層構造を持っているので、データ単位でランダムサンプリングして同一Sentinel Landscape内のデータを学習/テストに分けてCVしてしまうと、最終スコアリングに使われる学習/テストと異なる条件になってしまって、おそらく同一Sentinel Landscapes内でオーバーフィッティングしていたほうがいいスコアを出すデータセットができてしまう。なので、CVでは、同一Sentinel Landscape内のデータが学習とテストに分かれないように注意してfoldを分ける必要がある。
ちなみに各データがどのSentinel Landscapeに属するかを表すIDは振られていないので(データの説明に書くならふれよ!!と思った)、TMAPという項目をキーとして使う。TMAPは年間降水量かなにかの変数だけど、衛星から取ったデータの空間解像度の悪さからこの値をキーとしてデータをまとめるとほぼSentinel Landscapeの単位でまとめれるという議論がフォーラムであって、その話をしている人達はかなり信頼できる人達だったので、オレもこれを使うようにした。

コンペ後のフォーラムを見ると、上位のほうがみんなSentinel LandscapeベースのCVをしていたってことはないようだけど、ただ少なくともlocationベース(同じ位置から取られたSubsoil/Topsoilのデータが学習/テストで別れないようにする)でやらないと、オーバーフィッティングしてたほうが高いスコアが出てしまうので、順位落とした人の多くは間違った条件のCVでハイパーパラメータを調節した結果、元になったAbhishekのコードよりも悪い結果になってしまったのでは、と思います。

最後にpublic LBが信頼できないだろうと思った理由。

  • 計算に使っているデータが少なすぎる(90件くらい)
  • CVに比べて異常によいスコア(public LB: 約0.40, CV: 約0.50、しかもこの超簡単なデータセットはprivate LBの計算から除外されるので、private LBはいっそう難しくなる)
  • CVと負の相関すらあるように見える!!
  • 信頼できそうな人達がみんな信頼できないと言っていた(重要)

追記

なんとかツリー系、入力が信号で説明変数間に強い相関があるのと、次元に対して学習データが少なすぎるのとで、変数選択がうまくいかなくて、決定木をベースとしたものは良い結果が出せなかったのでは、と思っているけど詳しくないので分かりません。
ただ100位以内の多くは、SVRGBMの結果を混ぜたものみたいです。