幾何分布と負の二項分布
幾何分布は確率pのベルヌーイ試行において、試行回数がx回で初めて成功するときの分布でした。
f(x)=pqx−1where q=1−p
一方、負の二項分布は幾何分布の一般化として知られています。k回目の成功を得るまでの失敗の回数をxとすれば、
f(x)=(k−1k+x−1)pkqx
であります。ここでk=1であれば
f(x)=pqx
となります。教科書1では幾何分布の方の確率変数が成功までの試行回数に対して、負の二項分布の確率変数が失敗の回数であることに注意しないといけません。幾何分布の確率変数を失敗の数だと思って、最初のx回は全部失敗したという風に再解釈すれば、互いに分布はk=1のときに等しくなります。
次のコードでk = 1のときの幾何分布と負の二項分布を比較してみましょう。
import kotlin.math.exp
import kotlin.math.pow
import kotlin.math.min
import kotlin.math.ln
fun factorial(n: Int, acc: Int = 1): Int {
if (n == 0) return acc
return factorial(n - 1, n * acc)
}
fun Int.logCombination(k: Int): Double {
if (k < 0 || k > this) return Double.NEGATIVE_INFINITY
val kk = min(k, this - k)
return (1..kk).sumOf { i ->
ln((this - i + 1).toDouble()) - ln(i.toDouble())
}
}
fun geometricFromFailures(p: Double, x: Int): Double {
val q = 1.0 - p
return p * q.pow(x.toDouble())
}
fun negativeBinomialDistribution(k: Int, p: Double, x: Int): Double {
val q = 1.0 - p
return exp((k + x - 1).logCombination(x)) * p.pow(k) * q.pow(x)
}
fun main() {
val k = 2
val p = 0.5
val n = 10
val geometricMean = (1.0 - p) / p
val geometricStd = kotlin.math.sqrt((1.0 - p) / (p * p))
val nbMean = k * geometricMean
val nbStd = kotlin.math.sqrt(k * (1.0 - p) / (p * p))
val xData = (0..n).toList()
val negativeBinomialData = xData.map{ negativeBinomialDistribution(k, p, it) }
val geometricData = xData.map { geometricFromFailures(p, it)}
geometricData
.zip(negativeBinomialData)
.mapIndexed { x, (g, nb) ->
Triple(x, g, nb)
}
.forEach { (x, g, nb) ->
println(
String.format(
"x=%2d:\tg=% .6f,\tnb=% .6f,\tdiff=% .6e",
x, g, nb, g - nb
)
)
}
}
また実際にプロットすると下記のようなイメージになります。
参考資料
リポジトリ