Ahlfors "Prove the formula of Gauss"

After you have established that the RHS and the LHS differ by a multiplicative constant, all you are left to do it plug in $z=1$. If you pair up the factors in the RHS as $$\Gamma \left( \frac{1+k}{n} \right) \leftrightarrow \Gamma \left( \frac{n-k-1}{n} \right) ,$$ and apply the reflection formula $\Gamma(z) \Gamma(1-z)=\frac{\pi}{\sin \pi z}$, things will be easier IMO.

Edit:

Once you apply the reflection formula, you will have to deal with a product of sines. Please see this question in order to handle it.

Another common approach is to derive it from the limit definition of the gamma function. (See below.)

The multiplication formula can be written in the form

$$n^{nz-1/2} \prod_{k=0}^{n-1} \Gamma \left(z+\frac{k}{n} \right) = (\sqrt{2 \pi})^{n-1} \Gamma(nz)$$

Using the limit definition of the gamma function, we have

$$ \Gamma \left(z +\frac{k}{n} \right) = \lim_{m \to \infty} \frac{m! \ m^{z+\frac{k}{n}-1}}{(z+\frac{k}{n})(z+\frac{k}{n}+1) \cdots (z+\frac{k}{n} + m -1)}$$

Then using Stirling's formula, we get

$$ \begin{align} \Gamma \left(z+\frac{k}{n} \right) &= \lim_{m \to \infty} \frac{\sqrt{2 \pi m} (\frac{m}{e})^m m^{z+\frac{k}{n}-1}}{(z+\frac{k}{n})(z+\frac{k}{n}+1) \cdots (z+\frac{k}{n} + m -1)} \\ &=\lim_{m \to \infty} \frac{\sqrt{2 \pi} (\frac{mn}{e})^m m^{z+\frac{k}{n}-1/2}}{(nz+k)(nz+k+n) \cdots (nz+k + mn -n)} \end{align}$$

So

$$ n^{nz-1/2} \prod_{k=0}^{n-1} \Gamma \left(z+\frac{k}{n} \right) $$

$$ = n^{nz-1/2}\lim_{m \to \infty}\frac{(\sqrt{2 \pi})^{n} (\frac{mn}{e})^{mn} m^{nz-n/2} m^{\frac{1}{n} \sum_{k=1}^{n-1} k}}{(nz)(nz+1)\cdots (nz+n-1)(nz+n) \cdots (nz+mn-n)\cdots(nz+mn-1)} $$

$$ = \lim_{m \to \infty} \frac{(\sqrt{2 \pi})^{n} (\frac{mn}{e})^{mn} (mn)^{nz-1/2}}{(nz)(nz+1)\cdots (nz+n-1)(nz+n) \cdots (nz+mn-n) \cdots(nz+mn-1)} $$

Replacing $mn$ with $m$ shouldn't change the value of the limit (I think).

Therefore,

$$ \begin{align} n^{nz-1/2} \prod_{k=0}^{n-1} \Gamma \left(z+\frac{k}{n} \right) &=\lim_{m \to \infty} \frac{(\sqrt{2 \pi})^{n} (\frac{m}{e})^{m} m^{nz-1/2}}{(nz)(nz+1)\cdots(nz+m-1)} \\ &=\lim_{m \to \infty} \frac{(\sqrt{2 \pi})^{n-1}m! \ m^{nz-1}}{(nz)(nz+1)\cdots(nz+m-1)} \\ &= (\sqrt{2\pi})^{n-1}\Gamma(nz) \end{align}$$

EDIT:

Wikipedia states that the limit definition is

$$ \Gamma(t) = \lim_{n \to \infty} \frac{n! \ n^{t}}{t(t+1) \cdots (t+n)}$$

But notice that

$$ \Gamma(t-1) = \frac{\Gamma(t)}{t-1} = \lim_{n \to \infty} \frac{n! \ n^{t-1}}{(t-1)t \ldots (t+n-1)}$$

$$\implies \Gamma(t) = \lim_{n \to \infty} \frac{n! \ n^{t-1}}{t(t+1) \ldots (t+n-1)}$$