乱数あれこれゲームではあまり乱数を使ってはいけないと言われますが、そうは言っても乱数が無いと面白くないのもまたゲーム。そんな乱数についてのTIPSです。円内均一乱数ある円の中に点をランダムに穿ちたい事が結構あります。3.141592f;のようにするとこれは中心に近い座標ほど高確率で取得される不均一な円内乱数になってしまいます:こうではなくて、例えば半径1の円周に点が10個あるとした時に、半径2の円周に点が20個あると均一密度になりますよね。半径2の円周の方が半径1の円周よりも2倍選択されやすければ、後はその円周状のどこか1点を選ぶだけで円内均一乱数となります。半径に対する円周の選択されやすさをグラフ化するとこうなります:このままだと実はあまりうまく使えないので、グラフの下部分の面積が1になるように、[グラフの値]/[面積]という計算をしてあげます。すると上のグラフは、となります。見た目変わっていませんが、縦軸のスケールが倍になっていて、グラフの下部分の三角形の面積は1になっています。いう「面積が1のグラフ」は「確率密度分布」という名前が付いています。上のような三角形の形をしているものは「三角分布」と呼ばれます。bの範囲の面積」を求めるとその値が取られる確率が出てくると言う便利なものなんです。グラフの面積というとこれは積分のお話です。(0≦x≦1)と2次関数になる事がわかります。x=aとするとF(a)は0≦x≦aという範囲が取られる確率となります。上の三角分布の累積確率密度関数は次のような2次関数のグラフになります:さて、今多分私達が持っている生成可能な乱数は一様乱数です。縦軸は確率で、必ず0〜1に向かう右上がりのグラフとなります。この確率をランダムに選択し、それに対応する横軸の値を取ると、それがそのまま偏りのある乱数生成になってしまうんです:これは数学では「逆関数を求める」という作業になります。上のF(x)=x^2の逆関数Finv(x)は解析的に解く事ができて、となります:このグラフは横軸が0〜1の一様乱数、縦軸が取るべき円周の長さとなるわけです。ところで、円周の長さというのは半径に比例していて、2πrで求められるのでした。上のグラフの縦軸はそのまま半径と捉えても同じ事になります。えば、ある方向を中心軸として放射線状に何かをばらまきたい事があるとしましょう。その時、ばらまく角度が限定していて、かつその角度内で均等にばらまきたいねぇと思うわけです:こういう角度限定の球面均等乱数を考えてみます。球面の座標を表す方法はいくつかありますが、極座標で表すと色々つぶしが利きます。上のような片面がざっくり切れた球を真横から見るとこんな感じになるかと思います:あるθに対応したX軸上の点xを取る事ができて、後はその点xを中心としてYZ平面に平行な半径r(θ)の円の円周上のどこかの点を取ると球面上の1点を表す事ができます。cos(η)となりそうです。mはθが取る最大角度で0≦m≦πです。η(イータ)はx座標が決まった後の円の中心からの任意な角度です。yz座標を決めるのに使うだけです。一見するとこれで良さそうに見えるのですが、このままだと残念ながら不均一な分布になってしまいます。?例えばθが1度くらいの時を想像してみて下さい。ηがどんな値であっても随分とX軸近辺に点が近寄ってしまします。円の半径が小さいので点が密であるわけです。θが90度くらいの時は円の半径が大きい分点が疎になってしまいます。円の半径が大きい所ほどより多く選択されなければ均一な分布にならないんです。この感覚は�@の円内均一分布と似ていますよね。=120度を仮定したものです。横軸がθ(度)、縦軸がその時の半径r(θ)です。よって90度の時に一番多く選択されるような乱数を振れば球面一様になります。このグラフは確率分布になっているわけです。ここからは�@と同じ作業をしていきます。上のグラフを面積が1になる確率密度分布f(θ)に直します。その為には上の面積で縦軸の値を割り算すれば良いのでした。θが0〜m度の間での面積は定積分で表せるので、と計算できます。式としてはグラフの縦軸のスケールが1-cos(m)倍になっただけです。F(θ)が出てきます:これをグラフ化してみます:右肩上がりでθ=120度の所で確率1になっていますね。この累積確率密度関数の縦軸を一様乱数で取り、それに対応する横軸θを取ってあげると、半径r(θ)がうまい頻度で抽出される事になり、結果として球面一様乱数が得られると言うわけです。F(θ)の逆関数Finv(θ)は次のように解析的に解く事ができます:xは0〜1の一様乱数で、Finv(x)には0〜mまでの角度θが返ります。cos(η)θを求める所だけが補正されているわけです。実際うまくいってるのかなぁと思いDirectXで描画テストしてみたのがこちら:上行が補正前、下行が補正後で、左列がZ軸側から見た物で、右列が穴のあいているX軸側から見た物です。これを見ると補正前の分布は円半径が小さい所(右端)に点が集中しています。右列を見るとX軸を貫くあたりに密集しているのが分かると思います。角度θに対する半径r(θ)の確率密度関数等は式だけ示すと、となりまして、一番下段の逆関数がθを求める式となります。cos(η)となります。値が大きくなるにつれて頻度が比例して減り、ある最大値以上は頻度0になる以下のような分布を「三角分布」と言います:このグラフは0〜1の範囲で面積が丁度1になるので確率密度分布になっています。いう偏りを持った三角乱数を作ってみます。上のグラフから確率密度関数と累積確率密度関数は次のように導けます:累積確率密度関数F(x)の逆関数は次のように算出できます:最後の導きは皆大好き2次関数の解の公式です。2次関数なので解が2つ出てきますが、Finv(x)の範囲は0〜1なのでマイナス側のみ採用となります。rand()は0〜1までの一様乱数です。triRandomとなるわけです。円形集中分布ある点の回りに集中的に分布させたい事もあります。�Bの三角関数はx=0側程頻度が高いのである意味集中分布と言えます。2Dの円の中心から外に向かうにつれて疎になるような集中分布乱数を作るのは簡単で、実は�@でもう出てきています。rand()とすると自然と集中分布になります。その様子は�@で図として示してあります。半径rがより内側になるよう選択するには、一様乱数のべき乗を取るのが一つの手です。rand()とするとさらに円の中心点に集中するようになります。以下にn=2,4,6,8,10,12乗で全く同じ乱数列によってそれぞれ点を200点穿った図を掲載します:左上が2乗、右上が4乗、左真ん中が6乗と続きます。2乗と12乗(右下)を比べると同じ乱数列でも12乗の方がより中心点に寄っているのがわかります。0.5の時は丁度円形一様乱数になります。0.5にすると円周側により集中する面白い乱数にもなるわけです。最小値と最大値を持ち、その中央値が最大頻度となるような乱数が欲しい事が割とあります。正規分布乱数を使うと似たような事は出来るのですが、正規分布は出力される値が無限小〜無限大までなので最小値〜最大値という設定が出来ません。三角分布を使えばまぁ簡単に作れはするのですが(^-^;、ここでは2次関数を使ってチャレンジしてみます。以下のような2次関数な確率分布を想定します:最小値が0、最大値が1で最頻値(中央値、平均値)が0.5である確率分布です。上の2次関数は実は確率密度関数に既になっていまして、以下の式で表されます:ここから累積確率密度関数F(x)は、と算出できます。ここまでは簡単(^-^)。乱数化するにはF(x)の逆関数を求めなければなりません。いう3次方程式を解きます。3次方程式の解の公式はありますが、そのままだと「負数の累乗根」という頭がこんがらがる計算をしなくてはならなくなるため、別の作戦で攻めてみます。3次方程式はx^2の項があると面倒臭いので「カルダノの方法」で用いられる変換を行い、x^2項を無くします。下のyの3次方程式となります。ここでyの項をp、定数項をqとし、さらにp=-3a^2、q=-b*a^2と置き換えると、Wikipediaにある「ビエタの解」を用いる事が出来ます。ビエタの解は次の通りです:これで3次方程式の解が計算できるのですから驚きです。先の式から、aとbは、と計算できるので、ビエタの解に代入して一般式を得られます。nの値で解が3つ出てきますが、0≦y≦1に収まるのはn=2の時です:F(x)の所に0〜1の乱数を振ると、2次関数確率分布に従った0〜1までの乱数が出力されます。後は最小値と最大値にスケールしてあげれば良いだけです。実際上式の乱数をExcelでシミュレートすると次のような頻度分布が得られました:オレンジが理論値で青色のが10000個位サンプルした頻度分布です。大丈夫そうですね(^-^)。