Im trying to find the asymptotic expression for the beta function $B(x,y) = \frac{ \Gamma(x) \Gamma(y)}{\Gamma(x+y) } $. Using Stirling's approximation $\Gamma(x) \sim \sqrt{\frac{2\pi}{x} } (\frac{x}{e})^{x}$ it's easy to get that for $x,y \gg 1$, $$B(x,y) \sim \sqrt{2\pi} \frac{x^{x - \frac{1}{2}} y^{y - \frac{1}{2}}}{(x+y)^{x + y - \frac{1}{2}}} \equiv G(x,y).$$
However, we can also use Laplace's method: consider $$B(Ma,Mb) = \int_{0}^{1} z^{Ma} (1-z)^{Mb} dz = \int_{0}^{1} e^{Ma \times \ln(z) + Mb \times \ln(1-z)}dz,$$ the function $Ma \times \ln(z) + Mb \times \ln(1-z)$ has a unique maximum at $ z = \frac{a}{a+b}$ so by Laplace's method we can approximate the integral as $M \to \infty$ with:
\begin{align*}\int_{0}^{1} e^{Ma \times \ln(z) + Mb \times \ln(1-z)} dz & \sim \sqrt{\frac{2\pi}{M}} (\frac{a}{a+b})^{Ma}(\frac{b}{a+b})^{Mb} \frac{1}{\sqrt{\frac{(a+b)^{2}}{a} + \frac{(a+b)^{2}}{b}}} \\ & = \sqrt{\frac{2\pi}{M}} a^{Ma + \frac{1}{2}} b^{Mb + \frac{1}{2}} (a+b)^{-Ma -Mb - \frac{3}{2}}.\end{align*}
Replacing $Ma=x,Mb=y$ we get:
$$B(x,y) \sim \sqrt{2\pi} x^{x + \frac{1}{2}} y^{y + \frac{1}{2}} (x+y)^{-x-y - \frac{3}{2}} = \frac{xy}{(x+y)^{2}} G(x,y) = \frac{1}{2 + \frac{x}{y} + \frac{y}{x}} G(x,y).$$
So we suddenly get a coefficient that depends on the ratio of $x,y$ as they tend to infinity, whereas earlier we had no such coefficient. I've so far haven't been able to figure out where this difference is coming from since both methods should certainly yield the same result. I would really appreciate it if someone can help me figure out what I'm missing.