# How can I use the Stirling's approximation to approximate factorials?

I'd like to exploit Stirling's approximation during the symbolic manipulation of an expression. Essentially, I want replace Factorial[n] with n^n E^-n Sqrt[2 \Pi n] everywhere in a large expression.

How is it possible to achieve that?

As an example, I'd like to approximate a sum like this

$$\frac{n!}{(n+m_1)! (n-m_1)!} + \frac{n!}{(n+m_2)! (n-m_2)!} + \frac{n!}{ (n+m_3)! (n-m_3)!}$$ to $$\frac{F(n)}{F(n + m_1)Fun(n - m_1)} + \frac{F(n)}{F(n + m_2) F(n - m_2)} + \frac{F(n)}{F(n + m_3) F(n - m_3)},$$ where $$F(n) := \sqrt{2 \pi n} \, e^{-n} n^n$$ is the Stirling approximation of $$n!$$.

• Capital N is a built-in symbol, so better use n instead. Using ReplaceAll with the rule Factorial[n] -> n^n E^-n Sqrt[2 \[Pi] n] should do, no? Commented Feb 24, 2019 at 11:37
• Normal[Series[n!, {n, Infinity, 0}]]//Simplify does the job. This is written in the help to n!. Commented Feb 24, 2019 at 12:37
• Can you provide an example expression? All kinds of tricks available to manipulate combinatorial expressions. Commented Feb 24, 2019 at 13:20
• How about using /. Factorial[x_] -> x^x E^-x Sqrt[2 \[Pi] x]? Commented Feb 24, 2019 at 13:24
• It might actually be better to re-express everything in terms of Gamma[] (using FunctionExpand[] if necessary) and then use a replacement rule for Stirling. That way, you can handle Beta[], Binomial[], FactorialPower[], Pochhammer[]... Commented Feb 24, 2019 at 14:40

Note: updated code as per comments below. I believe this works:

Subscript[m, 1] = 2;
Subscript[m, 2] = 3;
Subscript[m, 3] = 4;
myf1[n_] =
n!/((n - Subscript[m, 1])! (n + Subscript[m, 1])!) +
n!/((n - Subscript[m, 2])! (n + Subscript[m, 2])!) +
n!/((n - Subscript[m, 3])! (n + Subscript[m, 3])!);
myf2[n_] = (myf1[n] /. (n_)! -> (Sqrt[2 Pi n] Exp[-n] n^n))
myf2[n] // TeXForm
myf1[100] // N
myf2[100] // N


$$\frac {e^n (n - 4)^{\frac {7} {2} - n} n^{n + \frac {1} {2}} (n + 4)^{-n - \frac {9} {2}}} {\sqrt {2 \pi }} + \frac {e^ n (n - 3)^{\frac {5} {2} - n} n^{n + \frac {1} {2}} (n + 3)^{-n - \frac {7} {2}}} {\sqrt {2 \pi }} + \frac {e^n (n - 2)^{\frac {3} {2} - n} n^{n + \frac {1} {2}} (n + 2)^{-n - \frac {5} {2}}} {\sqrt {2 \pi }}$$

2.923181964355016*10^-158
2.925623558480809*10^-158

• Indeed it works. I'm not familiar with the double underscore __. Why is this necessary here? Commented Sep 10, 2019 at 11:52
• The double (and triple) underscores mean take any expression (n, n+3, f(n), Sqrt[n], etc), and replace n in these expressions with the replacement value. I think anyway. Commented Sep 10, 2019 at 13:20
• No, double underscore means you can take/match a sequence of comma-separated items! You can use a single underscore everywhere in the code above. See reference.wolfram.com/language/ref/BlankSequence.html Commented Sep 11, 2019 at 12:03
• Ok, thanks. I updated my code above and indicated the changes. Commented Sep 11, 2019 at 16:01

As mentioned in the comments, just use Series:

series = Series[
n!/(n+m1)!/(n-m1)!+n!/(n+m2)!/(n-m2)!+n!/(n+m3)!/(n-m3)!,
{n, Infinity, 2}
];
series //TeXForm


$$\exp \left((1-\log (n)) n+O\left(\left(\frac{1}{n}\right)^5\right)\right) \left(\frac{3 \sqrt{\frac{1}{n}}}{\sqrt{2 \pi }}+\frac{\left(-4 \sqrt{2} \operatorname{m1}^2-4 \sqrt{2} \operatorname{m2}^2-4 \sqrt{2} \operatorname{m3}^2-\sqrt{2}\right) \left(\frac{1}{n}\right)^{3/2}}{8 \sqrt{\pi }}+O\left(\left(\frac{1}{n}\right)^{5/2}\right)\right)$$

If you want a normal expression, just use Normal:

Normal @ series //TeXForm


$$e^{n (1-\log (n))} \left(\frac{-4 \sqrt{2} \operatorname{m1}^2-4 \sqrt{2} \operatorname{m2}^2-4 \sqrt{2} \operatorname{m3}^2-\sqrt{2}}{8 \sqrt{\pi } n^{3/2}}+\frac{3}{\sqrt{2 \pi } \sqrt{n}}\right)$$

• Thanks! I commented on this in the comments to the questions. Commented Mar 5, 2019 at 14:52
• This answer is very interesting, but I believe is incomplete: I appreciate that Series[n!, {n, Infinity, 0}] // Normal gives the correct answer for $n!$. Now I'd say that $n!/(n-m)!/(m+n)!$ is approximated by Fun[n]/Fun[n + m1]/Fun[n - m1], with Fun[n_] := E^-n (1/n)^(-(1/2) - n) Sqrt[2 \[Pi]]. But this is different than Normal[Series[n!/(n + m1)!/(n - m1)!, {n, Infinity, 2}]] as this answer suggests, isn't it? What am I missing? Thanks. Commented Mar 22, 2019 at 18:18

Another way of doing what you ask, illustrating the use of Block. Define the expression in which to do the replacements

expr = Sum[n!/((n + m[i])! (n - m[i])!), {i, 1, 3}]
(* n!/((n - m[1])! (n + m[1])!) + n!/((n - m[2])! (n + m[2])!) +
n!/((n - m[3])! (n + m[3])!) *)


Do the replacement

Block[{Factorial = Sqrt[2 π #] Exp[-#] #^# &}, expr]
(* (E^n n^(
1/2 + n) (n - m[1])^(-(1/2) - n + m[1]) (n + m[1])^(-(1/2) - n -
m[1]))/Sqrt[2 π] + (
E^n n^(1/2 + n) (n - m[2])^(-(1/2) - n + m[2]) (n + m[2])^(-(1/2) -
n - m[2]))/Sqrt[2 π] + (
E^n n^(1/2 + n) (n - m[3])^(-(1/2) - n + m[3]) (n + m[3])^(-(1/2) -
n - m[3]))/Sqrt[2 π] *)

• For clarification: And this would work because n! is interpreted as Factorial[n] making the replacement become (Sqrt[2 Pi #] Exp[-#] #^# &) [n]? (Not sure if parenthesis necessary or not though) Commented Sep 14, 2019 at 18:44