How to understand $\Gamma(\alpha,1)$ distribution when $\alpha$ is large? It is said for large $\alpha$, $\operatorname{Gamma}(\alpha,1) \approx \operatorname{Normal}(\alpha,\alpha)$. Can anyone tell me why?
We know the mgf of a $\operatorname{Gamma}(\alpha, \beta)$ distribution is given by
$$M_X(t)=\left(\frac{1}{1-\beta t}\right)^\alpha, \quad t<\frac{1}{\beta}.$$
Hence, given $\beta=1$, we have mgf $\left(\frac{1}{1- t}\right)^\alpha$, $t<1$.
Then $$M_X{\left(\frac{t-\alpha}{\sqrt{a}}\right)}=\left(\frac{1}{1- \frac{t-\alpha}{\sqrt{a}}}\right)^\alpha, \quad t<1.$$
How to get the limit when $\alpha\to\infty$ ?
 A: Heuristically speaking, since a $\Gamma(n,1)$ is the sum of $n$ independent exponential variables with mean $1,$ and since the exponential distribution with mean $1$ has variance $1$,  the central limit theorem says that this will be approximately normal for large $n$, with mean $n$ and variance $n.$ More precisely, $ (S_n-n)/\sqrt{n}$ will converge in distribution to a standard normal, where $S_n$ is the sum of $n$ independent exponentials and thus has a $\Gamma(n,1)$ distribution.
Somewhat more formally, the $\Gamma(\alpha,1)$ distribution has pdf $f_X(x) = \frac{1}{\Gamma(\alpha)}x^{\alpha-1}e^{-x},$ from which we can compute a characteristic function $\phi_X(t) = E(e^{itX})=(1-it)^{-\alpha}.$ Then, shifting and rescaling, the characteristic function of $Z = (X-\alpha)/\sqrt{\alpha}$ is $$\phi_Z(t) = E(e^{it(X-\alpha)/\sqrt{\alpha}}) = e^{-it\sqrt{\alpha}}(1-it/\sqrt{\alpha})^{-\alpha},$$ which in the limit of $\alpha\to\infty$ is $e^{-\frac{1}{2}t^2},$ exactly the characteristic function of a standard normal.

It's maybe a little bit more burdensome, but you don't need to go through characteristic functions... just transforming the PDF gives $$ f_Z(z) =\frac{\alpha^{\alpha-1/2}}{\Gamma(\alpha)}(1+z/\sqrt\alpha)^{\alpha-1}e^{-\alpha} e^{-z\sqrt \alpha}$$ for $z > -\sqrt\alpha.$ And then with the help of the Stirling formula $\Gamma(\alpha)\sim \sqrt{2\pi}\alpha^{\alpha-1/2} e^{-\alpha}$ you can compute a limit of $\frac{1}{\sqrt{2\pi}}e^{-z^2/2}.$

(Edit in response to request in the comments.) As for how to compute the limit in the MGF case, we want the limit of $e^{-t\sqrt\alpha}(1-t/\sqrt\alpha)^{-\alpha}.$ This is probably easiest seen by taking logs $$ \log(M) = -t\sqrt\alpha -\alpha\log(1-t/\sqrt\alpha).$$ The taylor expansion of $\log(1-x)$ at $x=0$ is $-x-x^2/2+O(x^3),$ so we get $\lim\log(M) = t^2/2,$ so $\lim M = e^{t^2/2}.$
Note in your edit you are not correctly handling the MDF of a shifted and rescaled variable, which is not the same thing as plugging in a shifted and rescaled $t.$ If $Z=(X-\mu)/\sigma,$ $$M_Z(t) = E(e^{tZ}) = E(e^{t(X-\mu)/\sigma}) = E(e^{Xt/\sigma}e^{-\mu t/\sigma})=e^{-\mu t/\sigma}E(e^{X(t/\sigma)})= e^{-\mu t/\sigma} M_X(t/\sigma).$$
