[R] Product of 1 - probabilities

[Gabor Grothendieck]

There are several arbitrary precision packages available: gmp (an
interface to the GNU multi-precision library on CRAN) and
bc (an R interface to the bc arbitrary precision calculator):
http://r-bc.googlecode.com

There are also packages providing R interfaces to two computer
algebra systems and they both support not only arbitrary
precision but also exact calculation:

http://rsympy.googlecode.com
http://ryacas.googlecode.com

> library(bc)
> 1 - (1-10^-75)^10
[1] 0
> bc("1 - (1-10^-75)^10")
[1] ".0000000000000000000000000000000000000000000000000000000000000000000000000100000000000000000000000000"


On Thu, May 21, 2009 at 10:15 AM, Mark Bilton <mcbilton at hotmail.com> wrote:
>
> I am having a slight problem with probabilities.
>
> To calculate the final probability of an event p(F), we can take the product of the chance that each independent event that makes p(F) will NOT occur.
> So...
> p(F) = 1- ( (1-p(A)) * (1-p(B)) * (1-p(C))...(1-p(x)) )
>
> If the chance of an event within the product occurring remains the same, we can therefore raise this probability to a power of the number of times that event occurs.
> e.g. rolling a dice p(A) = 1/6 of getting a 1...
> p(F) = 1 - (1- (1/6))^z
> p(F) = 1 - (1-p(A))^z tells us the probabiltity of rolling a 1 'at least once' in z number of rolls.
>
> So then to R...
>
> if p(A) = 0.01; z = 4; p(F) = 0.039
>
> obviously p(F) > p(A)
>
> however the problem arises when we use very small numbers e.g. p(B) = 1 * 10^-30
> R understands this value
> However when you have 1-p(B) you get something very close to 1 as you expect...but R seems to think it is 1.
> So when I calculate p(F) = 1 - (1-p(B))^z = 1 to the power anything equals 1 so p(F) = 0 and not just close to zero BUT zero.
> It doesn't matter therefore if z = 1*10^1000, the answer is still zero !!
>
> Obviously therefore now p(F) < p(B)
>
> Is there any solution to my problem, e.g.
> - is it a problem with the sum (-) ? ie could I change the number of bits the number understands (however it seems strange that it can hold it as a value close to 0 but not close to 1)
> -or should I maybe use a function to calculate the exact answer ?
> -or something else
>
> Any help greatly appreciated
> Mark
> -
>
>
>
>
>
> _________________________________________________________________
>
>
>        [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>