べき乗則分布を生成する乱数生成器?

StackOverflow https://stackoverflow.com/questions/918736

  •  06-09-2019
  •  | 
  •  

質問

C++ コマンドライン Linux アプリ用のテストをいくつか書いています。べき乗則/ロングテール分布を使用して一連の整数を生成したいと考えています。つまり、いくつかの数字は非常に頻繁に取得されますが、ほとんどの数字は比較的まれです。

理想的には、 rand() または stdlib ランダム関数の 1 つで使用できる魔法の方程式があれば十分です。そうでない場合は、使いやすい C/C++ のチャンクが最適です。

ありがとう！

役に立ちましたか？

解決

これ Wolfram MathWorldのページ では、一様分布 (ほとんどの乱数発生器が提供するもの) からべき乗則分布を取得する方法について説明します。

短い答え（上記のリンクから導き出されました）：

x = [(x1^(n+1) - x0^(n+1))*y + x0^(n+1)]^(1/(n+1))

どこ y は一様変量であり、 n は分配力、 ×0 そして ×1 分布の範囲を定義し、 バツ はべき乗則分布変量です。

他のヒント

あなたは（確率分布関数（PDF）と呼ばれる）したい分布を知っていて、適切に正規化し、それを持っている場合は、

、あなたがして取得するCDF（可能な場合）を反転させ、累積分布関数（CDF）を得るためにそれを統合することができますあなたの希望に均一[0,1]分布から必要な転換ます。

だから、あなたが欲しい分布を定義することから始めます。

P = F(x)

（中xの[0,1]）を与えるために統合

C(y) = \int_0^y F(x) dx

これは反転させることができれば、あなたが得る

y = F^{-1}(C)

そこでrand()を呼び出し、最後の行にCように結果を差し込むとyを使用します。

この結果は、サンプリングの基本定理と呼ばれています。これは、正規化の要件および解析的機能を反転させる必要が面倒である。

別の方法として、あなたは拒否技術を使用することができますし、別の番号を投げるし、あなたの最初のスローによってindeicated場所でPDFと比較し、所望の範囲内に均一に番号を投げます。二スローはPDFを超えた場合に拒否します。長い尾を持つもののように、低確率領域の多いPDFファイルのために非効率的である傾向がある...

中間的なアプローチは、力ずくでCDFを反転含まれます。あなたは、ルックアップテーブルとしてCDFを格納し、その結果を得るために、逆引き参照を行う

。 <時間>

ここでの本当のstinkerは、単純なx^-n分布が範囲[0,1]上の非正規化可能であるということですので、あなたは、サンプリング定理を使用することはできません。試してみてください（X + 1）^ - Nの代わりに...

べき乗則分布を生成するために必要な数学についてはコメントできません (他の投稿に提案があります) が、TR1 C++ 標準ライブラリの乱数機能についてよく理解しておくことをお勧めします。 <random>. 。これらは、より多くの機能を提供します std::rand そして std::srand. 。新しいシステムは、ジェネレーター、エンジン、ディストリビューション用のモジュラー API を指定し、多数のプリセットを提供します。

含まれている配布プリセットは次のとおりです。

  • uniform_int
  • bernoulli_distribution
  • geometric_distribution
  • poisson_distribution
  • binomial_distribution
  • uniform_real
  • exponential_distribution
  • normal_distribution
  • gamma_distribution

べき乗則分布を定義するときは、それを既存の発電機やエンジンに接続できる必要があります。本 C++ 標準ライブラリ拡張機能 ピート・ベッカー著には、次のような素晴らしい章があります。 <random>.

ここに記事があります 他の分布の作成方法について (Cauchy、カイ 2 乗、Student t、および Snedecor F の例を含む)

私は（当然）受け入れ答えを補完するものとして、実際のシミュレーションを行ってみたかったです。 Rであるが、コードは-pseudoコード（擬似）となるように単純です。

おそらくより受け入れ答えおよびその他、でのWolfram MathWorld式の間の1つの小さな差一般的な、方程式は（典型的には、αとして示される）のべき乗則指数のnが明示負の符号を有していないという事実です。そう選択したアルファ値は、典型的には2〜3を負にする必要があります。

x0とx1は分布の下限と上限を表します。

だからここにある：

x1 = 5           # Maximum value
x0 = 0.1         # It can't be zero; otherwise X^0^(neg) is 1/0.
alpha = -2.5     # It has to be negative.
y = runif(1e5)   # Number of samples
x = ((x1^(alpha+1) - x0^(alpha+1))*y + x0^(alpha+1))^(1/(alpha+1))
hist(x, prob = T, breaks=40, ylim=c(0,10), xlim=c(0,1.2), border=F, 
col="yellowgreen", main="Power law density")
lines(density(x), col="chocolate", lwd=1)
lines(density(x, adjust=2), lty="dotted", col="darkblue", lwd=2)

" ここに画像の説明を入力する

または対数スケールでプロット

h = hist(x, prob=T, breaks=40, plot=F)
     plot(h$count, log="xy", type='l', lwd=1, lend=2, 
     xlab="", ylab="", main="Density in logarithmic scale")

" ここに画像の説明を入力する

ここではデータの要約である。

> summary(x)
   Min.   1st Qu.  Median    Mean   3rd Qu.    Max. 
  0.1000  0.1208  0.1584    0.2590  0.2511   4.9388
ライセンス： CC-BY-SA と 帰属
所属していません StackOverflow
scroll top