標準偏差を求めるperl関数


標準偏差っていうのが何なのか、最近知りました。

標準偏差 - Wikipedia

こういう数式で表されます。

\sigma^2 = \frac{1}{n} \sum_{i=1}^{n}(x_i - \bar{x})^2

ていうか、僕こういう数式にすごく抵抗あるんですよねぇ。
プログラムの方が理解しやすいかも。

これを見て、そういえばBuzzurlのタグクラウドを作るのに、標準偏差を基に重み付けをしたなあと思ってソースを見たら汚かったので書き直してみた。

普通に書くとこう。

#分散を求める。標準偏差が欲しい場合は戻り値をsqrt()
#average関数は自明なので省略
sub variance {
    my $num = scalar(@_) or return 0;
    my $ave = average(@_);
    my $ret = 0;
    while(@_){
        my $i = shift;
        $ret += ($i-$ave)**2;
    }
    return $ret/$num;
}

まあこれでいいんだけど長ったらしい。多少遅くなっていいとすれば、さらに短く書くほうがperlらしいかも？

#map/reduceを使うバージョン。
#ループが2回走る分、巨大な配列に対しては遅い。
use List::Util qw(reduce);
sub variance_mr {
    my $num = scalar(@_) or return 0;
    my $ave = average(@_);
    my $ret = reduce { $a + $b } map { ($_ - $ave) ** 2 } @_;
    return $ret/$num;
}

もっと遅くなってもいいなら、さらに短くこんなんでもいいかも。

#あらかじめmap内で濃度で割ってしまう版。
#割り算タイミングの問題で値が変わる場合がある。
sub variance_mr2 {
    my $num = scalar(@_) or return 0;
    my $ave = average(@_);
    return reduce { $a + $b } map { ($_ - $ave) ** 2 / $num } @_;
}

後者2つの書き方だと、元の数式とほとんど等価に見える気がする (( reduce { $a + $b } が \sum_{i=1}^{n}に相当する。
map { ($_ - $ave) ** 2 } @_が(x_i - \bar{x})^2に相当する。))
。あと、後者2つは遅いけど、MapReduceフレームワークみたいなものが使える人ならそれにのせやすいのでクソ分散実行できるかも？
まあ何が言いたいかというと、tex記法を使ってみたかっただけなんですが。