Coefficients of the stirling's series expansion for the factorial. Knowing the Stirling's approximation for the Gamma function (factorial) for integers:
$$\Gamma(n+1)=n!\approx \sqrt{2\pi n}n^ne^{-n}\bigg(1+\frac{a_1}{n}+\frac{a_2}{n^2}+\cdots\bigg)$$
Using the above approximation one can write:
$$(n+1)!=\sqrt{2\pi(n+1)}(n+1)^{n+1}e^{-(n+1)}\bigg(1+\frac{a_1}{n+1}+\frac{a_2}{(n+1)^2}+\cdots\bigg)$$
We know that following recursion holds:
$$(n+1)!=(n+1)n!$$
One can rewrite this:
$$(n+1)!=(n+1)\sqrt{2\pi n}n^ne^{-n}\bigg(1+\frac{a_1}{n}+\frac{a_2}{n^2}+\cdots\bigg)$$
All this comes from: https://www.csie.ntu.edu.tw/~b89089/link/gammaFunction.pdf (Page 8-9)
Then the author gives this expansion to calculate the $a_k$ coefficients when $n$ becomes large.
Comparing these two expressions for $(n+1)!$ gives
$$1+\frac{a_1}{n}+\frac{a_2}{n^2}+\cdots=\bigg(1+\frac{1}{n}\bigg)^{n+1/2}e^{-1}\bigg(1+\frac{a_1}{n+1}+\frac{a_2}{(n+1)^2}+\cdots\bigg)$$
Then he says, that after "classical series expansion" this equals:
$$1+\frac{a_1}{n}+\frac{a_2-a_1+\frac{1}{12}}{n^2}+\frac{\frac{13}{12}a_1-2a_2+a_3+\frac{1}{12}}{n^3}+\cdots$$
I don't understand how he got there. Only thing that came to my mind was, that as $n\to\infty$ $\big(1+\frac{1}{n}\big)^n$ goes to $e$ but then we are left with $\big(1+\frac{1}{n}\big)^{1/2}$. When i expand this into binomial series, i get:
$$\bigg(1+\frac{1}{2n}-\frac{1}{8n^2}+\frac{1}{16n^3}\cdots\bigg)\cdot\bigg(1+\frac{a_1}{n+1}+\frac{a_2}{(n+1)^2}+\cdots\bigg)$$
And I'm stuck here. 
Or is there any other elementary way how to compute the coefficients $a_k$ of the stirling's series expansion for factorial/Gamma function?

You were close. Note that1
$$
e^{-1}\left(1+\frac{1}{n}\right)^{n+1/2}=1+\frac{1}{12n^2}-\frac{1}{12n^3}+\frac{113}{1440n^4}+\cdots\tag1
$$
Therefore, you have
$$
\bigg(1+\frac{1}{12n^2}-\frac{1}{12n^3}+\cdots\bigg)\cdot\bigg(1+\frac{a_1}{n+1}+\frac{a_2}{(n+1)^2}+\frac{a_3}{(n+1)^3}+\cdots\bigg)=\\
1+\frac{a_1}{n}+\frac{a_2-a_1+\frac{1}{12}}{n^2}+\frac{a_3-2a_2+\frac{13}{12}a_1+\frac{1}{12}}{n^3}+\cdots\tag2
$$
as per
$$
\frac{1}{n+1}=\frac{1}{n}-\frac{1}{n^2}+\cdots\\
\frac{1}{(n+1)^2}=\frac{1}{n^2}-\frac{2}{n^3}+\cdots\\
\frac{1}{(n+1)^3}=\frac{1}{n^3}-\frac{3}{n^4}+\cdots\\
\text{etc.}\tag3
$$

Easy (but slightly cumbersome) exercise:
\begin{equation}
\begin{aligned}
{}&\lim_{n\to\infty}e^{-1}\left(1+\frac{1}{n}\right)^{n+1/2}=1\\
{}&\lim_{n\to\infty}n^2\left[e^{-1}\left(1+\frac{1}{n}\right)^{n+1/2}-1\right]=\frac{1}{12}\\
{}&\lim_{n\to\infty}n^3\left[e^{-1}\left(1+\frac{1}{n}\right)^{n+1/2}-1-\frac{1}{12n^2}\right]=-\frac{1}{12}\\
{}&\text{etc.}
\end{aligned}\tag4
\end{equation}
These are basic limits that are left to the reader.
