How to come up with the gamma function? It always puzzles me, how the Gamma function's inventor came up with its definition
$$\Gamma(x+1)=\int_0^1(-\ln t)^x\;\mathrm dt=\int_0^\infty t^xe^{-t}\;\mathrm dt$$
Is there a nice derivation of this generalization of the factorial?

In Leonhard Euler's Integral: A Historical Profile of the Gamma Function: In Memoriam: Milton Abramowitz by Philip J. Davis in The American Mathematical Monthly , Dec., 1959: Apparently, Euler, experimenting
with infinite products of numbers, chanced to notice that if n is a positive integer,
$$\small n!=\Bigg [ \left(\frac 2 1 \right)^n \frac1{n+1}\Bigg ]\,\Bigg [ \left(\frac 3 2 \right)^n \frac2{n+2}\Bigg ]\,\Bigg [ \left(\frac 4 3 \right)^n \frac3{n+3}\Bigg ]\cdots$$
and succeeding in transforming this infinity product into an integral, extending the factorial beyond integers, upon noticing that for certain values the infinite product yielded $\pi, $ suggesting areas of a circle.
But there is a really neat intuition already expressed in one of the answers, and beautifully presented by Robert Andrew Martin here, simply expanding the integral part in the integration by parts of a polynomial modulated by an exponential. In essence, the counterpart of the factorials in Taylor series.
For instance for $x^4$ (leaving constant of integration out):
$$\begin{align}\small
\int x^4\; e^{-x} \; dx &\small= -x^4\;e^{-x} +\int 4\,x^3\;e^{-x}\;dx\\
&\small = -x^4 \; e^{-x}-4\,x^3\;e^{-x} +\int_0^\infty 4\cdot 3\,x^2\;e^{-x}\;dx\\
&\small = -x^4 \; e^{-x} -4\,x^3\;e^{-x}  -4\cdot 3\,x^2\;e^{-x} +\int 4\cdot 3\cdot 2\,x\;e^{-x}\;dx\\
&\small =-x^4 \; e^{-x} -4\,x^3\;e^{-x}  -4\cdot 3\,x^2\;e^{-x}  -4\cdot 3\cdot2\,x\;e^{-x} - \underbrace{4\cdot 3\cdot 2\cdot 1}_{4!}\;e^{-x}
\end{align}$$
Generalizing and integrating from $0$ to $\infty:$
$$\small\int_0^\infty x^n\; e^{-x} \; dx=-x^n \; e^{-x} -n\,x^{n-1}\;e^{-x}  \cdots - \underbrace{n\cdot (n-1)\cdots 3\cdot 2\cdot 1}_{n!}\;e^{-x}\;\;\Bigr|_{x=0}^\infty=n!$$
which can immediately be extended beyond integers as $\displaystyle \small x! = \int_0^\infty t^x\; e^{-t} \; dt$
essentially the gamma function, except for the accepted slightly different definition: $$\Gamma(x)=\int_0^\infty t^{x-1}\; e^{-t} \; dt$$
that makes $\small (x-1)!=\Gamma(x).$
