# Stirling's approx for large factorials - need trick please

1. Sep 28, 2003

### mmwave

To calculate the multiplicities of 600 heads in 1000 coin tosses you start with 1000 choose 600 or
1000! / (600! * (1000-600)!) which equals 1000! / (600!) * (400!).

Since you can't calculate this easily, apply Stirling's approx.
N! = N^N * e^(-N) * sqrt(2piN). Applying this to numerator and denominator and leaving out sqrt terms since they are not the problem:

1000^1000 * e^(-1000)
--------------------------------
600^600*e-600 * 400^400 *e-400

Now the exponential terms cancel, but you still can't calculate the remaining terms.

Question: What is the trick to reduce this to something I can calculate? Here's my best attempt.

I tried factoring the bases and canceling terms but that gives

10*10 ^1000 * 10^1000 10^1000
--------------------------------------- = --------------------
10*10 ^600 * 10*10 ^400 * 6^400 * 4^400 6^600 * 4 ^400

You can cancel common factor of 2^1000 the same way giving

5^1000 / (3^600 * 2^400)

I still can't calculate this. Is there another approach or am I missing some other opportunity for canceling?

2. Sep 28, 2003

### Hurkyl

Staff Emeritus
If only there was a way to transform exponentiation into multiplication, and multiplication into addition...

3. Sep 28, 2003

### mmwave

I tried the log formula for large factorials first but I got numbers I couldn't take the antilog of. But you made me think of taking the logs where I left off and subtracting the log of 2^N then antilog. If that doesn't work then I'm still lost.

4. Sep 28, 2003

### Hurkyl

Staff Emeritus
Well, the result does have 291 decimal digits.

5. Oct 4, 2010

### mcglothlinj

ok.... you know omega=1000!/600!400!.
So ln(omega)=ln(1000!/600!400!)
sterling: N! \approx (N/e)^N*sqrt(2PiN)
so,
ln(omega)=ln((1000^1000)*sqrt(2000Pi)/((600^600)*sqrt(1200Pi)*(400^400)*sqrt(800Pi))

log rules:ln(ab)=lna+lnb , ln(a/b)=lna-lnb , ln(a^b)=blna

so,
ln(omega)=ln(1000^1000)+ln(sqrt(.02575161)-ln(600^600)-ln(400^400)
ln(omega)=1000ln1000-3.65925-600ln600-400ln400
ln(omega)=669.4
so, omega=e^669.4

let e^a=2^b, so ln(e^a)=ln(2^b), so a=bln2, so b=a/ln2

so omega=e^669.4=2^965.8
this is the number of ways to obtain 600 heads out of 1000 coins.

we know there exists 2^1000 possible outcomes from 1000 coin flips, so the probability is given by

P=omega/2^1000=2^965.8/2^1000=2^-34.2=5*10^-11

6. Oct 5, 2010

### awkward

The mean of the distribution of the number of heads is 500 and the standard deviation is $$\sqrt{600 * 0.5 * 0.5} \approx 15.8$$, so another approach is to evaluate the Normal pdf with that mean and standard deviation at x=600. This approximation gives p = 5.2 * 10^-11. By way of comparison, Wolfram Alpha says p = 4.6 * 10^-11 for 600 heads and 400 tails.

7. Oct 5, 2010

### CRGreathouse

Yeah, 16 microseconds with Pari, what a drag.

http://pari.math.u-bordeaux.fr/