Evaluating $\int_0^\infty x^{2\alpha-1}(e^{-x^\alpha}-e^{-r^\alpha x^\alpha})^{n-2}e^{-x^\alpha(1+r^\alpha)}\,dx$ I am looking for a tractable closed form of the following integral where $\alpha>0,r>1$:
$$I=\int_0^\infty x^{2\alpha-1}(e^{-x^\alpha}-e^{-r^\alpha x^\alpha})^{n-2}e^{-x^\alpha(1+r^\alpha)}\,dx\tag{1}$$
It can be assumed that $n\ge2$ is an integer.
Plugging $I$ directly into Mathematica, I get a sort of ugly answer involving non-elementary functions. Not to mention there is some confusion over the convergence of the integral.
Is there a way to tackle the integral by hand and arrive at a form involving elementary functions?
The context is in finding the density of $(\max_{1\le i\le n}X_i/\min_{1\le i\le n} X_i)$ where $X_1,X_2,\cdots,X_n$ are independent random variables having the common Weibull density 
$$f(x)=\alpha x^{\alpha-1}e^{-x^{\alpha}}\mathbf1_{x>0}\,,\quad\alpha>0$$
I had used a change of variables $(\max X_i,\min X_i)\to\left(\frac{\max X_i}{\min X_i},\min X_i\right)$ to be faced with $(1)$ trying to find the marginal density $\frac{\max X_i}{\min X_i}$ from the joint density of $\left(\frac{\max X_i}{\min X_i},\min X_i\right)$.
 A: You can start with the substitution $x^\alpha = y$ . If we let $R = r^\alpha$ , it leads to
$$ I = \frac{1}{\alpha} \int \limits_0^\infty y \left(\mathrm{e}^{-y} - \mathrm{e}^{-R y}\right)^{n-2} \mathrm{e}^{-(1+R) y} \,  \mathrm{d} y \, .$$
In this form it is even more evident that there are no convergence issues for the given ranges of the parameters.
If you want to avoid special functions, you should probably expand the integrand using the binomial theorem:
\begin{align}
I &= \frac{1}{\alpha} \sum \limits_{k=0}^{n-2} {n-2\choose k} (-1)^k \int \limits_0^\infty y \mathrm{e}^{-[R(k+1) + n-2-k+1] y} \,  \mathrm{d} y \\
&= \frac{1}{\alpha} \sum \limits_{k=0}^{n-2} {n-2\choose k} \frac{(-1)^k}{[n + (r^\alpha-1)(k+1)]^2} \, .
\end{align}
I do not see further simplifications for general values of the parameters, but it is a finite sum of nice functions.
In order to get a solution in terms of special functions, the substitution $\mathrm{e}^{-(R-1) y} = t$ can be used to find
\begin{align}
I &= - \frac{1}{\alpha(R-1)^2} \int \limits_0^1 \ln(t) t^\frac{n}{R-1} (1-t)^{n-2} \, \mathrm{d} t \\
&= - \frac{1}{\alpha (r^\alpha -1)^2} \operatorname{\partial_1} \operatorname{B} \left(\frac{n}{r^\alpha - 1} + 1 , n-1\right)
\end{align}
with the beta function $\operatorname{B}$, which eventually reduces to the result given by Mathematica.
