[R]疑似乱数生成(正規分布)

キーポイント確率統計 (理工系数学のキーポイント)」という本の「ポイント8確率過程とシミュレーション」 式(8.33, 8.34)(P。137-138) で説明されている 疑似正規乱数の生成を R で実装。実際に手を動かしたほうが理解も深まるしね、、、

単純 バージョン(8.33)

sum(runif(12)) - 6
★実行結果★

> replicate(10, sum(runif(12)) - 6)
 [1]  1.6246508 -0.7546210  0.7789311 -0.2014089  1.0188684
 [6]  0.6849198  0.1321929 -0.9739276 -1.3647944  0.1087528

Box-Muller Transform バージョン

box_muller  replicate(10, box_muller())
 [1]  0.6989254  1.7411153 -1.6949648 -0.4014906 -2.2848479
 [6] -2.2337043 -1.6000773  1.0380264  0.6702963 -0.2575914

ヒストグラムで分布を確認

まとめて実行
乱数の生成からグラフの作成まで。

#!/usr/bin/Rscript

print ('built-in normal distribution')
for (i in 1:10) {
  print (rnorm(1))
}

print ('simple normal distribution generation')
for (i in 1:10) {
  print (sum(runif(12)) - 6)
}

print ('box-muller transform generation')
box_muller <- function() {
  sqrt(-2 * log (runif(1))) * cos( 2 * pi * runif(1))
}
for (i in 1:10) {
  print (box_muller())
}

# make graph:START
png(file="normal-distribution.png")
par(mfrow=c(3, 1))
hist(rnorm(10000),                          main='built-in normal distribution',   br=100)
hist(replicate(10000, sum(runif(12)) - 6 ), main='simple normal distribution',     br=100)
hist(replicate(10000, box_muller()),        main='box-muller normal distribution', br=100)
dev.off()
# make graph:END

★実行結果★

$ Rscript  pseudo-normal-distribution.R
[1] "built-in normal distribution"
[1] -1.368874
[1] 0.9321853
[1] 0.6318225
[1] -0.7848614
[1] -0.8923767
[1] 1.804753
[1] -0.5966181
[1] 1.850231
[1] -0.3963262
[1] -0.7763592
[1] "simple normal distribution generation"
[1] 0.02792411
[1] 1.653229
[1] -1.493716
[1] -1.472915
[1] -0.04638198
[1] 1.629666
[1] -0.3880242
[1] -1.677303
[1] 0.03766184
[1] -0.3427567
[1] "box-muller transform generation"
[1] 0.7826891
[1] -1.417955
[1] -0.8144024
[1] 0.2500512
[1] 1.709031
[1] -0.1271732
[1] 0.5425739
[1] 0.1091272
[1] -0.6947928
[1] 1.055268
null device
          1
Advertisements
Tagged with: , , , ,
Posted in algorithm, R

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

Archives
  • RT @__apf__: How to write a research paper: a guide for software engineers & practitioners. docs.google.com/presentation/d… /cc @inwyrd 4 months ago
  • RT @HayatoChiba: 昔、自然と対話しながら数学に打ち込んだら何かを悟れるのではと思いたち、専門書1つだけ持ってパワースポットで名高い奈良の山奥に1週間籠ったことがある。しかし泊まった民宿にドカベンが全巻揃っていたため、水島新司と対話しただけで1週間過ぎた。 それ… 5 months ago
  • RT @googlecloud: Ever wonder what underwater fiber optic internet cables look like? Look no further than this deep dive w/ @NatAndLo: https… 5 months ago
  • @ijin UTC+01:00 な時間帯で生活しています、、、 10 months ago
  • RT @mattcutts: Google's world-class Site Reliability Engineering team wrote a new book: amazon.com/Site-Reliabili… It's about managing produc… 1 year ago
%d bloggers like this: