How to calculate the probability of N normal variable being in increasing order?

Question

Suppose we have $n$ normal variable $X_1,X_2,\dots,X_n$, with corresponding mean $\mu_1,\dots,\mu_n$ and sd $\sigma_1,\dots,\sigma_n$. What is the probability of $X_1 < X_2 < \dots < X_n$, i.e. $P(X_1 < X_2 < \dots < X_n)$.

Numerical method is also okay. Actually I have thought about how to do the multi-variable integration but I cannot get it done. When $n=2$, it is easy to see that $X_1-X_2$ is also a normal variable so $P(X_1-X_2<0)$ is quite easy to calculate. However, I think when $n=3$, because of $P(X_1-X_2)$ and $P(X_2-X_3)$ are not independent, the result is not easy to calculate. Currently I am thinking about whether we can use some method to approximately calculate the value the that probability. I need to come up with an algorithm to compute that.

Thank you!

Answer 1

We can use multivariate normal distribution to solve this problem. Denote the difference between consecutive elements with a (n-1)-dimensional random vector $Y=(X_1-X_2,X_2-X_3,\dots, X_{n-1}-X_n)$, what we need to find is $Y<(0,0,...,0)$.

Random vector $Y$ has a multivariate normal distribution.

The mean vector $\mu = (E(X_1-X_2), \dots, E(X_{n-1}-X_n)) = (\mu_1 - \mu_2, \dots, \mu_{n-1} - \mu_n)$.

The covariance matrix $\Sigma = [Cov[X_i-X_{i+1}, X_j-X_{j+1}]], i=1,\dots, n-1, j=1,\dots, n-1$.

For each entry of the matrix,

$Cov_{i,i} = Var(X_i-X_{i+1}) = \sigma_i^2 + \sigma_{i+1}^2$

$Cov_{i,i+1} = Cov_{i+1,i} = E((X_i-X_{i+1})(X_{i+1}-X_{i+2})) - E(X_i-X_{i+1})E(X_{i+1}-X_{i+2})$

$ = E(X_iX_{i+1} + X_{i+1}X_{i+2}-X_{i+1}^2-X_iX_{i+2})-(\mu_i-\mu_{i+1})(\mu_{i+1}-\mu_{i+2})$

Since $X_i, X_{i+1}, X_{i+2}$ are independent, we have the above expression

$= E(X_i)E(X_{i+1})+E(X_{i+1})E(X_{i+2})-E(X_{i+1}^2)-E(X_i)E(X_{i+2}) - (\mu_i-\mu_{i+1})(\mu_{i+1}-\mu_{i+2})$

Apart from $E(X_{i+1}^2)$, other terms are easy to calculate now. For $E(X_{i+1}^2)$, it can be calculate using chi-square distribution of order one. The result turns out to be $\mu_{i+1}^2+\sigma_{i+1}^2$. Plug in the above equation, we have $Cov_{i,i+1}=-\sigma_{i+1}^2$.

Other entries in the covariance matrix is zero.

Finally, use numerical method to calculate the probability. One implementation of the cdf of multivariate normal distribution can be found in http://www.math.wsu.edu/faculty/genz/software/software.html, the MVNDST function.