文系プログラマが数式アレルギーを治す為に、中学〜高校レベルの知識でもわかる数式を、プログラムのコードに直す作業を黙々と続ける企画。
ソースコードはJuliaを利用。5年後にはこの分野の主流言語になっているのではないかと勝手に思っているJuliaを利用。
今回のテーマは標準偏差。私の数学の偏差値は53くらいです。たぶん。
式はWikipediaの標準偏差とかを参考にした。
http://ja.wikipedia.org/wiki/%E6%A8%99%E6%BA%96%E5%81%8F%E5%B7%AE
こんな感じ。
# 自分用メモ $ \sigma^2 = \frac{1}{n} \sum\limits_{i=1}^n ({x_i - \mu})^2 $
上の式自体は分散を求める式だけど、分散の平方根を取れば標準偏差になるので気にしない。
σはシグマと読む。標準偏差を意味する変数に使われることが多い。
μ は ミューと読む。平均値を表現する時に使われることが多い。
Σ はシグマと読む。大文字シグマ。σは小文字シグマ。Σは総和を表す。上の例だと、iは1〜nまでの数列で、Σの部分はx(1) - μの2乗 〜 x(n) - μの2乗までの総和を意味する。
Juliaのコードに直すと、下記のようになる。平均値(μ)はmeanでさくっと取ることにする。
x = [3.0, 2.5, 1.4, 6.2, 3.3, 2.1, 5.0] μ = mean( x ) var_sum = 0 for i = x var_sum += ( i - μ ) ^ 2 end n = length(x) # 分散 var_sum / n #=> 2.4367346938775514 # 標準偏差 sqrt( var_sum / n ) #=> 1.5610043862454555
Juliaでは変数にアルファベットや数値以外も使えるので、μをそのまま変数名として利用できる。この分野で力を発揮する箇所だと思う。
上の処理では、まずはmeanで平均値(μ)を取って、あとは配列の中身を回しながら、(i - μ)の2乗の総和を出して、最後に配列のlengthで割っている。
ただlengthで割るだけだと分散に、sqrtで平方根を取ると標準偏差になる。
実際にJuliaの分散と標準偏差を出す機能を使って答え合わせをしておく。
x = [3.0, 2.5, 1.4, 6.2, 3.3, 2.1, 5.0] var( x ) #=> 2.8428571428571434 std( x ) #=> 1.6860774427223513
あれ、さっき出した結果と違うじゃん。
これはアレやね、この手の機能では標本に対する標準偏差を出す用であって、その場合は式がちょっと変わる。具体的にはnで割ってるところが、n-1で割ることになる。
# 自分用メモ $ \sigma^2 = \frac{1}{n-1} \sum\limits_{i=1}^n ({x_i - \mu})^2 $
では、実際にさっきのコードで、nで割っていたところを n - 1 として再度実行してみる。
x = [3.0, 2.5, 1.4, 6.2, 3.3, 2.1, 5.0] μ = mean( x ) var_sum = 0 for i = x var_sum += ( i - μ ) ^ 2 end n = length( x ) # 分散 var_sum / (n - 1) #=> 2.8428571428571434 var( x ) #=> 2.8428571428571434 # 標準偏差 sqrt( var_sum / (n - 1) ) #=> 1.6860774427223513 std ( x ) #=> 1.6860774427223513
よし、揃った。