日々之迷歩

世の中わからんことだらけ

ITが複雑で難しくなっていく様に翻弄される日々です。微力ながら共著させていただいた「シェル・ワンライナー160本ノック」をよろしくお願い申し上げます。

擬似正規乱数と確率密度分布

gnuplotが面白いので、調子に乗ってブログ2つ目。擬似正規乱数と確率密度分布のお話。

擬似一様乱数から擬似正規乱数を作成する方法として、12個ずつ足して6引くといいらしい。これもシェル芸本にやり方が載っていた。そもそもなぜこうすると擬似的な正規乱数になるのかは・・・別途統計学勉強・・・

gihyo.jp

データを扱う上でいろいろと統計学力が足りないことが分かっただけでも良しとし、気を取り直して擬似正規乱数の作成と、正規分布に似ているのか?確率密度分布のグラフを出してみることにした。

実行環境は下記の通り。

  • OSX El Capitan
  • GNU coreutils
  • GNU awk
  • ImageMagick
  • gnuplot
  • Open usp Tukubai

作ったシェルスクリプトはこちら。Nとdeltaを調節して、確率密度関数がどう変わるか実験出来る。実行するとdistribution.pngという画像ファイルに、確率密度関数のグラフが記載されている。

$ cat random.bash 
#!/bin/bash

N=1000000
delta=0.1
hdelta=$(echo $delta | awk '{print $1/2}')
area=$(echo $N $delta | awk '{print $1*$2}')

yes | gawk 'BEGIN{srand()}{print rand()}' | awk '{a+=$1;if(NR%12==0){print a-6;a=0}}' | ghead -n $N > uniform_random_number.txt

cat uniform_random_number.txt | awk "{i=0;if(\$1>0){i=int(\$1/$delta)*$delta}else{i=int(\$1/$delta)*$delta-$delta};print i}" | LANG=C gsort -n | count 1 1 | awk "{print \$1+$hdelta,\$2/$area}" > distribution.txt

cat << EOF > plot.gp
set terminal png size 800,800
set output 'distribution.png'
set nokey
set xr[-5:5]
set yr[0:0.5]
set title "Probability density distribution(Uniform Ramdom Number)\nN=$N: delta=$delta"; plot 'distribution.txt' with boxes
EOF

gnuplot plot.gp

$ ./random.bash
$ ls *.png
distribution.png

N=1000、N=10000、N=100000で実行した結果。Nが大きいほど正規分布に近づくことがわかる。これくらいのことならば、シェルスクリプトとgnuplotで十分なのかも。

f:id:papiro:20151220220653p:plain f:id:papiro:20151220220700p:plain f:id:papiro:20151220220707p:plain