On 5/21/2009 10:15 AM, Mark Bilton 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 ?
The problem is that R maintains about 15-16 decimal places of accuracy,
and to that accuracy, 1- 1.e-30 is equal to 1. You need to find a
different formula, e.g.
(1-x)^z = exp(z * log(1-x)) = exp(z * log1p(-x))
So if you have z = 1.e20, and p(A) = 1.e-30, then
1 - (1-1.e-30)^1.e20 can be calculated in R as
1 - exp(1.e20 * log1p(-1.e-30))
which works out to 1.e-10.
Duncan Murdoch
-or something else
Any help greatly appreciated
Mark
-
_________________________________________________________________
[[alternative HTML version deleted]]
______________________________________________
R-help@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.
______________________________________________
R-help@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.
