0
$\begingroup$

In my problem, $p$, $X$, and $N$ are fixed and I want to know the probability. Currently, I am calculating this recursively, as for any given $p$, $X$, $N$, the probability is:

function Prob(p,X,N) {
  return p*Prob(p,X-1,N-1) + (1-p)*Prob(p,X,N-1)
}

The intuition behind this is that if we were successful (first term) we have $X-1$ more successes to find and $N-1$ more trials. If we weren't successful $(1-p)$, we have $X$ more successes to find and $N-1$ more trials. There are also base cases that I left out where $X=0$ (probability is 1 because we've gotten all the successes) and $N$ less than $X$ (probability is zero because we have no more trials)
The recursion is really killing the speed of computation, so I'm looking for a closed form or approximation or just a point in the right direction.
Thanks

$\endgroup$
9
  • $\begingroup$ It is just the binomial term $_N C$$_X$ $p^X$ ($1-p$$^N$$^-$$^X$).. $\endgroup$ Feb 23, 2019 at 21:56
  • $\begingroup$ ah, sorry, it wasn't clear that it's x or more successes. that binomial equation gives me the equation to calculate the probability of exactly X successes $\endgroup$
    – Grifball
    Feb 23, 2019 at 22:03
  • 1
    $\begingroup$ Then just sum up the terms from X to N. $\endgroup$ Feb 23, 2019 at 22:05
  • $\begingroup$ that is a good idea, im trying that right now $\endgroup$
    – Grifball
    Feb 23, 2019 at 22:06
  • 1
    $\begingroup$ Thanks guys for helping out, it's running a lot faster after I implemented the 1 - cdf idea. $\endgroup$
    – Grifball
    Feb 23, 2019 at 23:44

1 Answer 1

0
$\begingroup$

Most programming languages have "Cumulative Distribution Functions" (CDF) for various distributions. As a result, you shouldn't have to individually compute each probability and then sum up to $X$. For example, in R, to obtain the CDF of a binomial distribution with $p=.20$, $n=50$, and $X=10$, you don't need to calculate:

> sum(dbinom(x=0:10, size=50, prob=0.20))
[1] 0.5835594

Instead, you could simply use:

> pbinom(q=10, size=50, prob=0.20)
[1] 0.5835594

Comparing the times, you see that the built in CDF function is generally faster:

> Sys.time()->start;
> sum(dbinom(0:54, size=100, prob=.50) )
[1] 0.8158992
> print(Sys.time()-start);
Time difference of 0.002012014 secs
> 
> Sys.time()->start;
> pbinom(q=54, size=100, prob=0.50)
[1] 0.8158992
> print(Sys.time()-start);
Time difference of 0.0009989738 secs

Also, if $n$ is large and/or $p$ is close to 0.50, you can use a normal approximation to the binomial. For example, suppose you wanted to calculate the Probability of at least 41 successes in 80 tries when $p=0.5$. You could calculate:

n=80
x=41
p=.5
q=1-p

pbinom(n,size=n,prob=p)-pbinom(x,size=n,prob=p)
[1] 0.3687772

or you could use the normal approximation:

1-pnorm(x,mean=n*p,sd=sqrt(n*p*q))
0.4115316

and with continuity correction, you get:

1-pnorm(x+.5,mean=n*p,sd=sqrt(n*p*q))
0.3686578
$\endgroup$
2
  • 1
    $\begingroup$ Wow! That is much faster, thanks. I'm using python, so the code ended up being: 1 - scipy.stats.binom.cdf(X-1, Y, p) I think my p's are going to very a lot so this should work for me. $\endgroup$
    – Grifball
    Feb 23, 2019 at 23:39
  • $\begingroup$ Great! Yes, most languages have CDF functions that are optimized for this so you don't have to "roll your own," which tend to be slower, in my experience. $\endgroup$ Feb 23, 2019 at 23:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.