# Approximate order statistics for normal random variables

Are there well known formulas for the order statistics of certain random distributions? Particularly the first and last order statistics of a normal random variable, but a more general answer would also be appreciated.

Edit: To clarify, I am looking for approximating formulas that can be more-or-less explicitly evaluated, not the exact integral expression.

For example, I have seen the following two approximations for the first order statistic (ie the minimum) of a normal rv:

$e_{1:n} \geq \mu - \frac{n-1}{\sqrt{2n-1}}\sigma$

and

$e_{1:n} \approx \mu + \Phi^{-1} \left( \frac{1}{n+1} \right)\sigma$

The first of these, for $n=200$, gives approximately $e_{1:200} \geq \mu - 10\sigma$ which seems like a wildly loose bound.

The second gives $e_{1:200} \approx \mu - 2.58\sigma$ whereas a quick Monte Carlo gives $e_{1:200} \approx \mu - 2.75\sigma$, so it's not a bad approximation but not great either, and more importantly I don't have any intuition about where it comes from.

Any help?

-
If you use R, see the ppoints function. –  cardinal Mar 31 '11 at 13:41
@probabilityislogic has given some good intuition for the approximations you list. Would it be helpful at all if I gave some more from an alternative viewpoint, or have you satisfied your curiosity on this matter? –  cardinal Apr 1 '11 at 19:59

The classic reference is Royston (1982)[1] which has algorithms going beyond explicit formulas. It also quotes a well-known formula by Blom (1958): $E(r:n) \approx \mu + \Phi^{-1}(\frac{r-\alpha}{n-2\alpha+1})\sigma$ with $\alpha=0.375$. This formula gives a multiplier of -2.73 for $n=200, r=1$.

[1]: Algorithm AS 177: Expected Normal Order Statistics (Exact and Approximate) J. P. Royston. Journal of the Royal Statistical Society. Series C (Applied Statistics) Vol. 31, No. 2 (1982), pp. 161-165

-
Thanks, this is exactly what I was looking for. –  Chris Taylor Mar 31 '11 at 21:25


1. $i-1$ values below $x$, this has probability $F_{X}(x)$ for each observation, where $F_X(x)=\Pr(X<x)$ is the CDF of the random variable X.
2. $N-i$ values above $x$, this has probability $1-F_{X}(x)$
3. 1 value inside a infinitesimal interval containing $x$, this has probability $f_{X}(x)dx$ where $f_{X}(x)dx=dF_{X}(x)=\Pr(x<X<x+dx)$ is the PDF of the random variable $X$

There are ${N \choose 1}{N-1 \choose i-1}$ ways to make this choice, so we have:

$$f_{i}(x_{i})=\frac{N!}{(i-1)!(N-i)!}f_{X}(x_{i})\left[1-F_{X}(x_{i})\right]^{N-i}\left[F_{X}(x_{i})\right]^{i-1}dx$$

EDIT in my original post, I made a very poor attempt at going further from this point, and the comments below reflect this. I have sought to rectify this below

If we take the mean value of this pdf we get:

$$E(X_{i})=\int_{-\infty}^{\infty} x_{i}f_{i}(x_{i})dx_{i}$$

And in this integral, we make to change of variables to $p_{i}=F_{X}(x_{i})$ (taking @henry's hint), and the integral becomes:

$$E(X_{i})=\int_{0}^{1} F_{X}^{-1}(p_{i})\Beta(p_{i}|i,N-i+1)dp_{i}=E_{\Beta(p_{i}|i,N-i+1)}\left[F_{X}^{-1}(p_{i})\right]$$

So this is the expected value of the inverse CDF, which can be well approximated using the delta method to give:

$$E_{\Beta(p_{i}|i,N-i+1)}\left[F_{X}^{-1}(p_{i})\right]\approx F_{X}^{-1}\left[E_{\Beta(p_{i}|i,N-i+1)}\right]=F_{X}^{-1}\left[\frac{i}{N+1}\right]$$

To make a better approximation, we can expand to 2nd order (prime denoting differentiation), and noting that the second derivative of an inverse is:

$$\frac{\partial^{2}}{\partial a^{2}}F_{X}^{-1}(a)=-\frac{F_{X}^{''}(F_{X}^{-1}(a))}{\left[F_{X}^{'}(F_{X}^{-1}(a))\right]^{3}}=-\frac{f_{X}^{'}(F_{X}^{-1}(a))}{\left[f_{X}(F_{X}^{-1}(a))\right]^{3}}$$

Let $\nu_{i}=F_{X}^{-1}\left[\frac{i}{N+1}\right]$. Then We have:

$$E_{\Beta(p_{i}|i,N-i+1)}\left[F_{X}^{-1}(p_{i})\right]\approx F_{X}^{-1}\left[\nu_{i}\right]-\frac{\Var_{\Beta(p_{i}|i,N-i+1)}\left[p_{i}\right]}{2}\frac{f_{X}^{'}(\nu_{i})}{\left[f_{X}(\nu_{i})\right]^{3}}$$ $$=\nu_{i}-\frac{\left(\frac{i}{N+1}\right)\left(1-\frac{i}{N+1}\right)}{2(N+2)}\frac{f_{X}^{'}(\nu_{i})}{\left[f_{X}(\nu_{i})\right]^{3}}$$

Now, specialising to normal case we have $$f_{X}(x)=\frac{1}{\sigma}\phi(\frac{x-\mu}{\sigma})\rightarrow f_{X}^{'}(x)=-\frac{x-\mu}{\sigma^{3}}\phi(\frac{x-\mu}{\sigma})=-\frac{x-\mu}{\sigma^{2}}f_{X}(x)$$ $$F_{X}(x)=\Phi(\frac{x-\mu}{\sigma})\implies F_{X}^{-1}(x)=\mu+\sigma\Phi^{-1}(x)$$

Note that $f_{X}(\nu_{i})=\frac{1}{\sigma}\phi\left[\Phi^{-1}\left(\frac{i}{N+1}\right)\right]$ And the expectation approximately becomes:

$$E[x_{i}]\approx \mu+\sigma\Phi^{-1}\left(\frac{i}{N+1}\right)+\frac{\left(\frac{i}{N+1}\right)\left(1-\frac{i}{N+1}\right)}{2(N+2)}\frac{\sigma\Phi^{-1}\left(\frac{i}{N+1}\right)}{\left[\phi\left[\Phi^{-1}\left(\frac{i}{N+1}\right)\right]\right]^{2}}$$

And finally:

$$E[x_{i}]\approx \mu+\sigma\Phi^{-1}\left(\frac{i}{N+1}\right)\left[1+\frac{\left(\frac{i}{N+1}\right)\left(1-\frac{i}{N+1}\right)}{2(N+2)\left[\phi\left[\Phi^{-1}\left(\frac{i}{N+1}\right)\right]\right]^{2}}\right]$$

Although as @whuber has noted, this will not be accurate in the tails. In fact I think it may be worse, because of the skewness of a beta with different parameters

-
"Maximum likelihood estimator of a random variable"? Not sure what that is, but I think you've (almost) calculated the mode. –  cardinal Mar 31 '11 at 14:39
Something mysterious happens about two-thirds of the way through when suddenly $\mu$ and $\sigma$ appear without warning or definition. –  whuber Mar 31 '11 at 15:17
I don't mean to "pile on", but it's also hard for me to see how the quantity in brackets can be approximated by a negative number. –  cardinal Mar 31 '11 at 15:34
mean(qnorm(qbeta(((1:100000) - 0.5 ) / 100000, 1, 200))) in R will give an approximate numerical integration, and in this case gives -2.746042. Just change the 1, 200 to get other quantiles. –  Henry Mar 31 '11 at 21:58
@cardinal - my mistake, this should be the inverse CDF $-\Phi^{-1}(\frac{1}{N+1})$ as defined in the question. And I have calculated the mode of the PDF, which is (but for an idealogical restriction) a maximum likelihood estimator - all I have to do is replace the "random variable" with a "greek letter" and then magically - its an MLE! But clearly my answer is in need of repair –  probabilityislogic Apr 1 '11 at 9:48

Depending on what you want to do, this answer may or may not help - I got the following exact formula from Maple's Statistics package.

with(Statistics):
X := OrderStatistic(Normal(0, 1), 1, n):
m := Mean(X):
m;


$$\int _{-\infty }^{\infty }\!1/2\,{\frac {{\it \_t0}\,n!\,\sqrt {2}{ {\rm e}^{-1/2\,{{\it \_t0}}^{2}}} \left( 1/2-1/2\, {{\rm erf}\left(1/2\,{\it \_t0}\,\sqrt {2}\right)} \right) ^{-1+n}}{ \left( -1+n \right) !\,\sqrt {\pi }}}{d{\it \_t0}}$$

By itself this isn't very useful (and it could probably be derived fairly easily by hand, since it's the minimum of $n$ random variables), but it does allow for quick and very accurate approximation for given values of $n$ - much more accurate than Monte Carlo:

evalf(eval(m, n = 200));
evalf[25](eval(m, n = 200));


gives -2.746042447 and -2.746042447451154492412344, respectively.

(Full disclosure - I maintain this package.)

-
@ProbabilityIsLogic derived this integral for all order statistics in the first half of his reply. –  whuber Mar 31 '11 at 19:42