2012年8月8日水曜日

一様分布から標準正規乱数を発生させる

標準正規分布に従う乱数を発生させる簡単な方法として、一様乱数をいくつか発生させてその和をとり、期待値で引き算するという方法がある。そこで、得られた値が本当に標準正規分布といって差し支えないか調べてみることにした。以下のコードでは、6つの一様乱数(取りうる値の範囲は[0,1])とし、それらの和の期待値3を引いたものを大量に発生させ、値の分布が一様乱数であるか検定したものである。6としているのは、とりうる値を[-3, 3]とするためである。標準正規分布であれば、多くの場合、そこから発生する乱数は[-3, 3]の範囲に収まるためである。

# 一様乱数から正規乱数を生成するプログラム

set.seed(1)
pts <- 6 # 一つの標準正規乱数を発生させるのに用いる一様乱数の数
N <- 5000 # 発生させる正規乱数
num.data <- N*pts
rn <- runif(num.data)
rnd.mat <- matrix(rn,nrow=N) # 一様乱数の行列
num <- rowSums(rnd.mat)-pts/2

hist(num,breaks=20,xlim=c(-4,4),xlab="",ylab="freq",main="")
oldpar <- par(new=T)
curve(dnorm,-4,4,xlab="",ylab="",axes=FALSE)
par <- oldpar

# コルモゴロフ・スミルノフ(KS)検定
# 標準正規分布 N(0, 1)とデータを比較している
ks.test(num, "pnorm", mean=0, sd=1)

コルモゴロフ・スミルノフ検定は、データがある分布に従っているかを検定するというものである。上のコードは1標本KS検定に相当する。帰無仮説はデータは分布から発生したものではない、対立仮説はデータは分布から発生したものとなる。統計量は従うと想定される分布の累積分布関数と、データのそれを比較した結果から得られる差の値に基づく。統計量が小さくなると、帰無仮説が棄却されるということになる。結論は、KS検定ではp値はほぼ0となり、帰無仮説は標準正規分布ではないということになった。ただし、中心極限定理からデータ数を6から増やせば、標準正規分布に従うことが想定される。そこで、pts<-12としてみると標準正規乱数といえるようなp値となった(有意水準5%で帰無仮説を棄却できる)。必要なptsのサイズは標準化する分布(今回は一様分布)の対称の程度に依存していることに注意。今回利用した一様分布U(0,1)は0.5で対称なので、比較的少ないデータ数で標準正規分布が生成できるのである。

0 件のコメント:

コメントを投稿

フォロワー

ページビューの合計