3

使用此代码时,我遇到了一些奇怪的问题:

positions<-c(58256)
occurrencies<-c(30)
frequency<-c(11/5531777)
length<-c(4)

prob<-c(0)
for(i in 0:(occurrencies-1))
{
  pow<-frequency^i
  pow1<-(1-frequency)^(positions-i)
  bin<-choose(positions, i)
  prob<<-prob+(bin*pow*pow1)
}

i此 for 循环的每次迭代都应计算在给定频率的情况下事件发生次数的二项式概率。每次迭代也会总结结果。这应该会导致prob变量永远不会超过 1,但是在 7 次左右的循环迭代之后,一切都变得糟糕并prob超过 1。

我认为这可能是数字精度的问题,所以我尝试使用Rmpfr但无济于事——同样的问题仍然存在。

我想知道是否有任何技巧或包来克服这种情况,或者我是否坚持这一点。

4

2 回答 2

3

你可以通过做来避免你的for循环

prob<-0
i    <- 0:(occurrencies-1)
pow  <- frequency^i
pow1 <- (1-frequency)^(positions-i)
bin  <- choose(positions, i)
prob <- cumsum(prob+(bin*pow*pow1))
[1] 0.8906152 0.9937867 0.9997624 0.9999932 0.9999998 1.0000000 1.0000000 1.0000000 1.0000000
[10] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
[19] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
[28] 1.0000000 1.0000000 1.0000000

我不知道这是否是您想要的结果,但您肯定可以避免这种for循环。

请参阅@Ben Bolker 的评论并查看pbinom功能。

于 2012-10-11T17:00:24.500 回答
2

按照 Ben Bolker 的建议,看看?pbinom

pbinom(q = occurencies, size = positions, prob = frequency, lower.tail = FALSE)
于 2012-10-12T17:24:34.157 回答