• ベストアンサー

三角分布に従う乱数

三角分布に従う乱数を発生させるプログラム 最頻値gで、0.8g~2.5gの範囲の三角分布(当然面積は1です)に従う乱数を発生させるプログラムを 書きたいと思っています。 この三角分布の確率密度関数をP(x)とすると、三角分布であるので最頻値gの左側である、傾きが正の直線h(x)と 最頻値gの右側である、傾きが負の直線f(x)で表せますよね 分布に従う乱数を発生させるためには、これら直線の関数を積分したものの逆関数x=P^-1(u) (Pのインバースです) (uは区間[0, 1]の一様乱数)とすればいいというところまでわかったんですが、 とりあえずh(x)とf(x)をそれぞれ積分して逆関数H^-1(u)、F^-1(u)を求めたところまではいいんですが x=H^-1(u)+F^-1(u)としてプログラムを実行すると、最頻値gの2倍あたりの値(例えば20に対して39など) しか出ず、最頻値gより小さい値が出ません。 H^-1(u)+F^-1(u)としているのがダメだと思うのですが、逆関数が2つある場合、ここからどうすればいいですか? また、初歩的な質問なのですが、区間[0, 1]の一様乱数というのはどう記述すればよいですか? ぜひ多くの方の回答お待ちしています。 よろしくお願いします。 (最頻値gは入力で与えるものとします)

質問者が選んだベストアンサー

  • ベストアンサー
noname#227064
noname#227064
回答No.2

最頻値がgで[a, b]の範囲で三角分布に従う乱数を発生させたいとします。 まず、確率密度関数を求めると、 (-∞, 0] : 0 (a, g] : 2(x-a)/{(g-a)(b-a)} (g, b] : 2(x-b)/{(g-b)(b-a)} (b, ∞) : 0 分布関数は確率密度関数を積分して、 (-∞, 0] : 0 (a, g] : (x-a)^2/{(g-a)(b-a)} (g, b] : (g-a)/(b-a)+(x^2-2bx-g^2+2bg)/{(g-b)(b-a)} (b, ∞) : 1 分布関数の逆関数を求めると、 (-∞, 0] : a (0, (g-a)/(b-a)] : a+sqrt((g-a)(b-a)x) ((g-a)/(b-a), 1] : b-sqrt((b-g)^2-(b-g)(b-a)(x-(g-a)/(b-a))) (1, ∞) : b となりますので、xは0以上1以下の値が入っているとすれば、 if (x <= (g-a)/(b-a)) return a+sqrt((g-a)(b-a)x); else return b-sqrt((b-g)^2-(b-g)(b-a)(x-(g-a)/(b-a))); とでも書けばいいでしょう。 (実際に書く場合は、もう少し計算量が少なくなるように書くとは思いますが) 計算があっているかどうかは確認しておいてください。

ou_gakusei
質問者

お礼

分かりやすい回答をどうもありがとうございました。 おかげでちゃんと実行することができました。 私のはそもそも積分が間違っていたようです。 詳しい解説、本当にありがとうございました。

その他の回答 (3)

  • Ishiwara
  • ベストアンサー率24% (462/1914)
回答No.4

ふつうの長方形分布の乱数を発生させ、条件に合わなければ捨ててしまう、というやり方が一番簡単です。

noname#227064
noname#227064
回答No.3

ANo.2 思いっきり間違いました。 こちらが正しいです。 if (x <= (g-a)/(b-a)) return a+sqrt((g-a)*(b-a)*x); else return b-sqrt((b-g)*(b-g)-(b-g)*(b-a)*(x-(g-a)/(b-a)));

noname#227064
noname#227064
回答No.1

言語は何を使われているのかわかりませんが、if文等で条件分岐をするのは駄目なのでしょうか? > また、初歩的な質問なのですが、区間[0, 1]の一様乱数というのはどう記述すればよいですか? rand関数というような関数がありませんか? 名前は違うかもしれませんが、必ず乱数を発生させる関数があるはずです。 ひょっとして、0からnまでの整数の乱数は得られるけど、[0, 1]の一様乱数をどうやって得るのかわからないということでしょうか? (C言語なら(double)rand() / RAND_MAXとでもすればいいですが)

ou_gakusei
質問者

お礼

回答ありがとうございます。 言語はC++です。乱数わかりました。 if文等で条件分岐をするというのは どういう条件のときにどういう動作をすればよいのかが難しいです・・

関連するQ&A

専門家に質問してみよう