Solution (source code)

= Solution

The given <probability density function> is that of a <beta distribution> with parameters $r$ and $n-r+1$. We can derive the required <moments> directly rather than assuming the beta moment formulas. For an integer $j\geq0$, the <Euler beta function> gives
$$
\begin{aligned}
\mathbb E[U_{(r)}^j]
&=\frac{n!}{(r-1)!(n-r)!}
\int_0^1u^{r+j-1}(1-u)^{n-r}\,du\\
&=\frac{n!(r+j-1)!}{(r-1)!(n+j)!}.
\end{aligned}
$$
Here $B(a,b)=(a-1)!(b-1)!/(a+b-1)!$ for positive integers $a,b$; this identity also follows by repeated <integration by parts>. Taking $j=1,2$ yields
$$
\mathbb E[U_{(r)}]=\frac r{n+1},\qquad
\mathbb E[U_{(r)}^2]=\frac{r(r+1)}{(n+1)(n+2)}.
$$
Subtract the square of the first <moment> from the second to obtain the <variance>:
$$
\boxed{\mathbb E[U_{(r)}]=\frac r{n+1},\qquad
\operatorname{Var}(U_{(r)})=
\frac{r(n-r+1)}{(n+1)^2(n+2)}.}
$$
These are exact finite-sample identities for each <uniform order statistic>.