「統計学関連なんでもあり」の過去ログ--- 040

No.03548　乱数の発生方法　　【初心者】　2007/06/03(Sun) 23:33

細菌121株につて最小阻止濃度（MIC：細菌の増殖を阻止する濃度）のデータがあります。

このデータに基づき，エクセルで乱数を発生させたいのですが，
正規分布に基づいた，NORMINV(RAND(),μ，σ)で発生させるのか，
対数正規分布に基づいて，＝EXP(NORMINV(RAND(),μ，σ）で発生させるのか？
どちらがよろしいでしょうか？
全くの素人で悩んでおります。よろしくお教え下さい。

No.03549　Re: 乱数の発生方法　　【青木繁伸】　2007/06/04(Mon) 08:30

どういう目的で乱数を発生させるのですか？
このデータと「同じ分布」をする乱数が欲しいのですか？
「全く同じ分布」をする乱数でも良いのですか？

No.03552　Re: 乱数の発生方法　　【初心者】　2007/06/04(Mon) 11:31

このデータと「同じ分布」をする乱数が欲いのですが，過去ログも参考にさせて頂きましたが，発生方法が分からなくて困ってます。何とぞお教えの程，お願いします。

No.03553　Re: 乱数の発生方法　　【青木繁伸】　2007/06/04(Mon) 11:59

このデータが対数正規分布に従うとして

MIC の，2を底とする対数をとる
株を重みとして，平均値μと分散σ^2（nで割る方）を計算する
2^NORMINV(RAND(),μ,σ) で元の尺度に戻す

なお，正規乱数を NORMINV(RAND(),μ,σ) で求めるのはお勧めできません

No.03556　Re: 乱数の発生方法　　【初心者】　2007/06/04(Mon) 16:38

ありがとうございます。今ひとつ理解できないのですが，

MICの対数（低2）をとり，μ＝-0.28926, σ=1.028251
株数(n=121)で重みつけして，

重み付μ=-0.28926/121=-0.00239
重み付σ=1.028251/121=0.008498

これを基に，
=EXP(NORMINV(RAND(), 重み付μ ,重み付σ)
で乱数を発生させるのでしょうか？

大変初歩的な質問で申し訳ありませんが，お教えの程お願いします。

No.03557　Re: 乱数の発生方法　　【青木繁伸】　2007/06/04(Mon) 17:24

わかりにくかったですか。
株を重みとしてというのは，たとえば 0.25 が4株ですから，LOG(0.25, 2) つまり-2 が 4個で，これを -2 の重みが 4 という意図でした。
（μ=1.190083, σ=1.654232 は，MIC の値をそのまま使ったものでしょう？これでは対数正規乱数を発生できません）
R なら
> MIC <- c(0.25, 0.5, 1, 2, 4, 8)
> n <- c(4, 54, 50, 5, 2, 6)
> log2MIC <- log(MIC, 2) # 対数変換
> ( mu <- sum(log2MIC*n)/sum(n) ) # 対数変換した後の平均値★
[1] -0.2892562
> ( sigma <- sqrt(sum((log2MIC-mu)^2*n)/sum(n)) ) # 対数変換した後の標準偏差●
[1] 1.023993
> result <- 2^rnorm(100000, mean=mu, sd=sigma) # 対数正規乱数
> # 平均値と，標準偏差が，元のデータを対数変換したものから計算したものと一致するか確認
> log2result <- log(result, 2) # 対数変換して
> mean(log2result) # 平均値 ★と近い値か
[1] -0.2948541
> sd(log2result) # 標準偏差 ●と近い値か
[1] 1.024733
なお，2を底とする必要はありません。自然対数でも，その他の底でもかまいません。

また，元のデータが，0.5 と 1 の所に極端に集中しているので，このデータから求めた平均値や標準偏差がどの程度の精度を持っているのかは相当に疑わしいでしょうね。

No.03558　Re: 乱数の発生方法　　【初心者】　2007/06/04(Mon) 17:45

青木先生　ありがとうございます。
何となくRの構文も少しですが見えてきました。

> result <- 2^rnorm(100000, mean=mu, sd=sigma) # 対数正規乱数

↑
ここで，resultに対数正規分布の値が格納されているのですね。
これをテキスト形式かCSV形式で書き出すことは出来ますでしょうか。

No.03559　Re: 乱数の発生方法　　【青木繁伸】　2007/06/04(Mon) 17:57

以下のあたりを参照されればよいかと
http://aoki2.si.gunma-u.ac.jp/R/io/index.html

一行に1つずつ書くなら
write(result, file="data", ncolumns=1)
などと