Stirling gave a formula for estimating n! but here is better |

Presuming an offset c improving the approximation, e ~ (n+1/n)^{(n+c)},
we estimate by reduction considering a numerically advanced
e ~ (n²/n²-1)^{(n²-1+c)},
pulling-out similar factors of e, = ((n/n-1)*(n/n+1))^{(n²-1+c)} ~
e^{(n+1-c)}*(n/n-1)^{(c²-c)} *
e^{(-n+c)}*(n/n+1)^{(c²+c-1)}
which has factor(s) that would be unity, and zeroing-out the
slowest, second term in (1-(1/n²))^{(c-c²)}*(1+(1/n))^{(1-2c)}
~ 1, whence c = 0.5, QED.

Alternatively, e ~ (n+.5/n-.5)^{n} for aesthetic symmetry. [Or faster,
(n+.5+(1/12n)/n-.5+(1/12n))^{n}]

(We could have used L'Hospital's Rule on log of numerator and denominator going to infinity: lim (c-c²)(log n-1) + (1-c-c²)(log n+1) / (1-2c²)(log n) = 1; by derivation, lim (c-c²/n-1) + (1-c-c²/n+1) / (1-2c²/n) = 1; by design faster to 1, whence -1 + 2c = 0.)

Stirling's approximation to n! ~ (n/e)^{n}*√2πn rewritten
(n/e)^{(n+.5)}*√2πe (which has a nice feature that there are no
stray factors of n) yields n+1!/n+1 ~
(n+1/e)^{(n+.5)}*√2πe^{-1} (which has a nice feature that it
has an estimate for 0!)... the exponent has not changed: only baseline parameters
n,√e have become n+1,√e^{-1}, exactly as in our approximation of
e (above). But midpoint peak (n+.5/e)^{(n+.5)}*√2π (*) is yet-nearer
~ n!, asserted as follows ...
(and also completely removes stray factors, and constants)

* (the slope of its curve, = 0 at c = 0:
δ_{c} (n+.5-c/e)^{(n+.5)}*√2π*e^{c} =
-(n+.5/e) + (n+.5-c/e) times a common positive: n ≠ -0.5)

[reconstruction ahead]

The comeback approximations to n = n!/n-1! are, by Stirling's rewritten
(n/e)(n/n-1)^{(n-.5)}, the next higher (n/e)(n+1/n)^{(n+.5)},
and midpoint
(n/e)√(2n+1/2n)^{(2n+.5)}(2n/2n-1)^{(2n-.5)}*(1-(1/4n²))^{.25},
closer by reason of the square-root being a geometric average of doubly-higher order
approximations to e; and the fourth-root being minuscule second-order compensation
~ (1/16n²) to even that (*). The initial 0! approximations are, Stirling's = 0,
next higher, √2π/e ~ 0.922, and midpoint, √(π/e) ~ 1.075, which starts
closest to unity,- and its successive factors, ×1×2×3×..., keep
it closest:-- this is a "sterling approximation" to n! .

* (In fact -(1/48n²) is even closer-still than -(1/16n²) but hors d'oeuvre.)

[A faster expansion takes the baseline (n+.5) to (n+.5-(1/24(n+.5))) -- note that (n+.5) then does triple duty]

This also yields an approximation to n! for noninteger (and fractional) n .

(We can more finely estimate the value of fractional factorials by n! = n+k!/(n+k)(n+k-1)...(n+1) : k → N>>0) .

By the approximation, δ n! ~ (n+.5/e)^{(n+.5)}*√2π*(log n+.5) ,
which comes back to n!*(log n+.5) ; And as (log n+.5)-(log .5) is the first
approximation to the series sum, Σ^{n} (1/i) ,
thus δ:0! ~ (log .5) = -(log 2) . Cf the gamma function.

A premise discovery under the title,

'Majestic Service in a Solar System'

Nuclear Emergency Management