I would like to write down an asymptotic expansion in the $N\to\infty$ limit of the following incomplete beta function
$$B\left(\frac{N}{N+1};N,p+1\right)=\int_0^{\frac{N}{N+1}}x^{N-1}(1-x)^p\,\text{d}x$$
where $N\in\mathbb{N}$ and $p>0$.
Clearly, when one deals with beta functions
$$B(a,b)\equiv B(1;a,b)=\frac{\Gamma(a)\Gamma(b)}{\Gamma(a+b)},$$
things are quite easy because one can expand the gamma functions by means of Stirling's approximation. Instead, I do not know how to proceed in my case due to the dependence on $N$ of the integration limit too.
While I was looking for a solution on the internet, I came across this article that gives the following asymptotic expansion for the incomplete beta function
\begin{equation}\tag{1}\label{eq1}
B(x;a,b)\sim\frac{x^a}{a}\sum_{k=0}^{\infty}\frac{f_k(b,x)}{a^k}\qquad\text{as}\;\;a\to\infty
\end{equation}
with $0\le x\le 1$ and
$$f_k(b,x)=\frac{\text{d}^k}{\text{d}w^k}\left[\left(1-x\text{e}^{-w}\right)^{b-1}\right]_{w=0}.$$
The formula is derived through a repeated integration by parts of the Laplace transform representation of the incomplete beta function
$$B(x;a,b)=x^a\int_0^{\infty}\text{e}^{-aw}\left(1-x\text{e}^{-w}\right)^{b-1}\text{d}w.$$
I decided to apply the above result to my specific case considering $N/(N+1)$ fixed, so after some simple algebra I obtained (truncating the series to $k=3$, I hope without errors)
\begin{equation}\begin{split}
B\left(\frac{N}{N+1};N,p+1\right)&\sim\frac{1}{N^{p+1}}\left[1+p+p(p-1)-\frac{p}{N}+p(p-1)(p-2)\right.\\[8pt]
&\left.\quad-\frac{3p(p-1)}{N}+\frac{p}{N^2}+\dots\right]\qquad\text{as}\;\;N\to\infty
\end{split}\end{equation}
First of all, I do not know if this approach is correct, especially considering that terms which are of order $O(1)$ in the square brakets appear for all $k$ of the above series.
I performed some simple simulations with Mathematica and the expansion seems to work well when $p\in\mathbb{N}$, but not for $p\not\in\mathbb{N}$. In particular, in this last case the truncation error taken from the article
$$\epsilon_K\equiv B(x;a,b)-\frac{x^a}{a}\sum_{m=0}^{K-1}\frac{f_m(b,x)}{a^m}$$
diverges as $K$ increases.
Edit: I just realized that the above expansion does not work for $p\not\in\mathbb{N}$ due to the derivation procedure when one integrates by parts the Laplace transform representation of the incomplete beta function. In fact, this produces a result's structure that allows for a spontaneous truncation of the asymptotic expression only for $p\in\mathbb{N}$.
I wonder if a result exists similar to the one of Eq. (\ref{eq1}) but valid for non-integer $b$.
Best Answer
Following a more elementary route we will first break the integral in two pieces:
$$B(\frac{N+1}{N+2};N+1,p+1)=B(N+1,p+1)-I\\I=\int^{1}_{\frac{N+1}{N+2}}x^{N}(1-x)^pdx$$
Perform the following change of variables $x=1-\frac{t}{N+2}$:
$$I=\frac{1}{(N+2)^{p+1}}\int_{0}^1t^p(1-\frac{t}{N+2})^{N}dt$$ However it is possible to estimate that $$(1-\frac{t}{N+2})^{N}=e^{-t}(1+\frac{4t-t^2}{2N}+\mathcal{O}(N^{-2}))$$
and performing the integrals we get:
$$I\sim\frac{1}{N^{p+1}}\Big[\gamma(p+1,1)+\frac{2}{N}(\gamma(p+2,1)-p-1-\frac{1}{4}\gamma(p+3,1))\Big]$$
While on the other hand we may write
$$B(N+1,p+1)\sim \frac{\Gamma(p+1)}{N^{p+1}}\Big[1-\frac{(p+1)^2+3(p+1)}{2N}\Big]$$
So all in all as a leading approximation one gets:
$$B(\frac{N+1}{N+2};N+1,p+1)\sim\frac{\Gamma(p+1)-\gamma(p+1,1)}{N^{p+1}}=\frac{\Gamma(p+1,1)}{N^{p+1}}+\mathcal{O}(N^{-p-2})$$
where $\Gamma(z,x)=\int_x^{\infty}t^{z-1}e^{-t}dt$ is the upper incomplete Gamma function. From the above one can calculate also the subleading term. More subleading terms can also be readily calculated but probably not in closed form.