D言語でコルモゴロフ-スミルノフ検定
この記事はD言語 Advent Calendar 2014 24日目の投稿です。
なぜこの題材を選んだのか
私は研究室内で行われたプログラミング講習にて、DDoS検知システムを作りました。その検知ロジックにコルモゴロフ-スミルノフ検定を利用したのですが、Rubyで実装したため、とても遅いのが問題でした。
そこで、ネイティヴで動き、C++よりもはるかに習得が容易だと私のまわりでHotなプログラミング言語であるD言語で実装し直してみることにしました。
コルモゴロフ-スミルノフ検定とは
コルモゴロフ-スミルノフ検定（以下、KS検定）とは、標本Xと標本Yの分布の相違を計る仮説検定のひとつです。
出典：http://sysplan.nams.kyushu-u.ac.jp/gen/edu/MarineStatistics/2007/kougi12/kougi12.pdf
計算手順
1. 標本Xの累積確率分布 
と、標本Yの累積確率分布
を求める。
2. 上記の2つの累積確率分布の差の絶対値の最大値であるKS統計量Dを求める。
3. Dから検定統計量
、相関係数pを算出する。
4. 帰無仮説の採否を決める(有意水準5%の場合)
p > 0.05の場合 : 帰無仮説の採択
p <= 0.05の場合 : 帰無仮説の棄却
カイ2乗分布を求める
qchi2はカイ2乗分布の上側累積確率を求める関数です。実装は以下の通りです。
//正規分布の下側累積確率
auto p_nor(immutable double z)
in {
assert(z >= 0.0 && z <= 3.09);
}
out(probability) {
assert(probability >= 0.0 && probability <= 0.499);
}
body {
double t, p;
t = p = z * exp(-0.5 * z * z) / sqrt(2 * PI);
for(size_t i = 3; i < 200; i+= 2) {
immutable auto prev = p;
t *= z * z / i;
p += t;
if (p == prev) { return 0.5 + p; }
}
return (z > 0.0) ? 1.0 : 0.0;
}
// 正規分布の上側累積確率
auto q_nor(immutable double z)
in {
assert(z >= 0.0 && z <= 3.09);
}
out (probability) {
assert(probability >= 0.001 && probability <= 0.5);
}
body {
return 1.0 - p_nor(z);
}
//上側累積確率
auto q_chi2(immutable int df, immutable double chi2)
in {
assert(df >= 1 && df <= 300);
assert(chi2 >= 0.00003927 && chi2 <= 366.844);
}
out(probability) {
assert(probability >= 0.0 && probability <= 1.0);
}
body {
double s, t;
if ((df & 1) == 1) { //自由度が奇数
const auto chi = sqrt(chi2);
if (df == 1) { return 2.0 * q_nor(chi); }
s = t = chi * exp(-0.5 * chi2) / sqrt(2 * PI);
for (size_t i = 3; i < df - 1; i += 2) {
t *= chi2 / i;
s += t;
}
return 2 * (q_nor(chi) + s);
} else { //自由度が偶数
s = t = exp(-0.5 * chi2);
for (size_t i = 2; i < df - 1; i += 2) {
t *= chi2 / i;
s += t;
}
return s;
}
}
KS検定を行う
KS検定を行うプログラムは以下の通りになります。
//累積和を求める
auto cumsum(immutable double[] a) {
auto b = cast(double[])a;
foreach(i; 1 .. a.length) {
b[i] = a[i] + a[i - 1];
}
return b;
}
auto ks_test(immutable double[] x, immutable double[] y)
in {
assert(x.length == y.length);
}
out(results) {
assert(results[2] >= 0.0 && results[2] <= 1.0);
}
body {
immutable auto sum_x = x.reduce!("a + b");
immutable auto sum_y = y.reduce!("a + b");
auto cum_x = cumsum(x).map!(s => s / sum_x);
auto cum_y = cumsum(y).map!(s => s / sum_y);
double[] cum = new double[cum_x.length];
for(size_t i = 0; i < cum.length; i++) {
cum[i] = abs(cum_x[i] - cum_y[i]);
}
// 累積分布関数の差の最大値
immutable auto d = cum.reduce!(max);
// 検定統計量
immutable auto chi = 4 * d * d * sum_x * sum_y / (sum_x + sum_y);
// 相関係数(p値)
immutable auto p = [1.0, pow(q_chi2(2, chi), 2.0)].reduce!(min);
return [d, chi, p];
}
void main() {
immutable auto x = [1.0,2.0,1.0,3.0,2.0,3.0,3.0,2.0,7.0,11.0,10.0,9.0,13.0,13.0,22.0,17.0,23.0,20.0,17.0,14.0,13.0,5.0,2.0,1.0,1.0,1.0];
immutable auto y = [0.0,1.0,2.0,2.0,4.0,5.0,6.0,8.0,10.0,7.0,16.0,17.0,17.0,13.0,19.0,13.0,18.0,10.0,4.0,6.0,6.0,5.0,1.0,3.0,0.0,0.0];
immutable auto p = ks_test(x, y)[2];
// 帰無仮説の採否(有意水準5%)
// p <= 0.05 (帰無仮説の棄却：二つの分布に差があると言える)
// p > 0.05 (帰無仮説の採択：二つの分布に差があると言えない)
(p <= 0.05).writeln();
}
実装は以上となります。全体の実装は以下のサイトに置いてあります。
感想
顔文字プログラミング言語を作るよ
この記事はKazoo04 Advent Calendar 2014 21日目の記事です。
どういうものなの？
顔文字でプログラミングができます。こちらで書かれたLisp処理系を使わせていただきました。ソースコードはgithubに公開してます。
https://github.com/olanleed/kaomoji_lang
作ろうと思ったきっかけ
確かかずー氏のとあるツイートがきっかけだったのだが、そのツイートがどこにあるのかわからない。そのツイートは私にとって特別だったのだろうか？ そもそもそんなツイートがあったのだろうか？ でもそれは思い出せないだけで、きっと私の心の奥底にある。
仕組み
処理系のコードを弄るのは手間がかかるのでやりたくはありません。しかし私は絶世の超絶天才爆殺エリートプログラマーなので、顔文字で書かれたソースコードをgsubメソッドでLispコードに置換することにより、秒速でこのプログラミング言語を創造することができました。まさに神技です。
言語仕様
まだ割り当ててないものがありますが、仕様は以下の通りです。
| Lisp | 顔文字言語 |
|---|---|
| () | ( ˘⊖˘) 。o() |
| + | 三(^o^)ノ |
| - | ▂▅▇█▓▒░(’ω’)░▒▓█▇▅▂ |
| * | ( ✹‿✹ )開眼だァーーーーーーーーーーー！！！！！！！！！！！ |
| / | ( (☛(◜◔。◔◝)☚)) |
| = | L(՞ਊ՞)｣ |
| modulo | (`o´) |
| not | (╮╯╭) |
| > | (:3[＿＿] |
| < | [＿＿]ε:) |
| >= | (¦3[＿＿] |
| <= | [＿＿]ε¦) |
| lambda | ( ◠‿◠ )☛ |
| display | (´へεへ`*) ＜ |
| define | (´へωへ`*) |
| if | +。:.ﾟ٩(๑＞◡＜๑)۶:.｡+ﾟ |
サンプルコード
(5 - 2) * 3を計算して結果を出力するプログラムです。
( ˘⊖˘) 。o((´へεへ`*) ＜( ˘⊖˘) 。o(( ✹‿✹ )開眼だァーーーーーーーーーーー！！！！！！！！！！！ ( ˘⊖˘) 。o(▂▅▇█▓▒░(’ω’)░▒▓█▇▅▂ 5 2) 3))
本当は階乗を求めるプログラムを書きたかったのですが、それさえしんどくなるほど記述力に優れた言語です。
Q & A
これは実用的な言語ですか？
実用的な言語の話をしている
gsubで置換していくとかパフォーマンス最悪では
これに何を求めてるんだ
今後の展望は
ないです
俺がkazoo04だ
質問ではない
Dropoutの実装と重みの正則化
この記事はMachine Learning Advent Calendar 2013 3日目の記事です.
はじめに
ニューラルネットワークの汎化性能を向上させるDropoutは, Deep Learningを実装する上で必須の技術だと思います. 本日はDropoutとその実装方法について説明させていただきます.
Dropoutとは
ニューラルネットは複雑なモデルであるため過学習に陥りやすいです. これを回避するためにはL2ノルムで値の増加を防いだり, L1ノルムでスパースにしたりするのが一般的です.
しかし正則化でもニューラルネットのような複雑なモデルに適切に制約を加えるのは困難です.
そこでDropoutの考え方です. Dropoutは各訓練データに対して中間素子をランダムに50%（または任意の確率で）無視しながら学習を進めます. 推定時は素子の出力を半分にします.
なぜこれだけで汎化性能が向上するのかというと, "アンサンブル効果"で説明がつきます.
アンサンブル効果
アンサンブル学習
アンサンブル学習とは個々に学習した識別器を複数用意し, 出力の平均を取るなど, それらをまとめあげて一つの識別器とする方法です. 一人で考えるよりも複数人で議論したほうが良い結果を得られるのと同じでしょうか.
Dropoutの実装
いろいろなアプローチがあると思いますが, ここではドロップアウトしたい隠れ層の素子をマスクする方法を紹介します.
Back Propagate
Back Propagate時も素子をなかったことにしないといけません.
もし活性化関数にシグモイド関数を使っていた場合は以下のように式を変形することができます.
これでDropoutの実装は終わりです！
重みの正則化
よりよい重みを速く求める, そして更に汎化性能を向上させるには重みの正則化をしたほうがよいとされています.
Momentum
モーメント法と呼ばれる勾配法の高速化手法です. 学習率が低いとなかなか収束しないため, 重みを調整してあげるのがmomentumの役割です.
勾配法によるBackPropagationの重み更新式は以下のように与えられています.
: 学習率
これにmomentum係数を導入して式を以下のように修正します.
: momentum係数
Weight decay
汎化性能向上のための手法の一つで, Weight decayを導入すると重みの爆発が抑制され, 過学習を防いでくれます.
: Weight decay係数
他にも中間素子に入力する重みのL2ノルムが閾値を越えたら, 閾値と同じ値となるように線形スケーリングを行う方法もあります[1].
Dropoutの欠点
Dropoutの欠点として, たくさん更新しないと学習が進みません（計算時間が増える）.
この欠点を克服した"Fast Dropout"がありますが, これはまた別の機会に紹介したいと思います.
おわりに
BackPropagationにDropoutを適用したソースコードを公開します. 間違いがありましたらご連絡ください.
