Obtaining a positive definite covariance matrix of order statistics

Suppose $X_1,\dots,X_n$ are independent samples from some distribution with known absolutely continuous CDF $F:\mathbb{R}\rightarrow[0,1]$. Let $X_{(1)},\dots,X_{(n)}$ denote the order statistics, i.e. the ordered sample. Defining the column vector $X_{(\cdot)}=[X_{(1)},\dots,X_{(n)}]'$, we want to numerically calculate $\mathbb{E}[(X_{(\cdot)}-\mathbb{E}X_{(\cdot)})(X_{(\cdot)}-\mathbb{E}X_{(\cdot)})']$.

The most obvious approach uses the fact that for $j,k\in\{1,\dots,n\}$, $j<k$:

$$\mathbb{E}X_{(k)}=\int{x\frac{n! [F(x)]^{k-1}[1-F(x)]^{n-k}}{(k-1)!(n-k)!} dF(x)},$$ $$\mathbb{E}X_{(k)}^2=\int{x^2\frac{n! [F(x)]^{k-1}[1-F(x)]^{n-k}}{(k-1)!(n-k)!} dF(x)},$$ $$\mathbb{E}X_{(j)}X_{(k)}=\int{\int{xy\frac{n! [F(x)]^{j-1}[F(y)-F(x)]^{k-1-j}[1-F(y)]^{n-k}}{(j-1)!(k-j-1)!(n-k)!} 1[x\le y] dF(x) } dF(y)},$$

where $1[\cdot]$ is the indicator function. See e.g. Wikipedia here for an informal proof.

Using these formulae with standard numerical integration methods works well for small $n$. However, for large $n$ the resulting covariance matrix often ends up non-positive definite unless implausibly many integration nodes are used. This is despite the individual elements of the covariance matrix usually being (loosely) close to a covariance matrix generated via a naïve Monte Carlo approach that guarantees positive definiteness (i.e. draw such a sample of length $n$, then sort it, repeat this lots of times, take the covariance).

Is it possible to express these integrals in such a way that the resulting covariance matrix is guaranteed to be positive definite?

E.g. is it the case that:

$$\mathbb{E}[(X_{(\cdot)}-\mathbb{E}X_{(\cdot)})(X_{(\cdot)}-\mathbb{E}X_{(\cdot)})']=\int{\int{g(u,v) dF(u)}dF(v)},$$

for some function $g:\mathbb{R}^2\rightarrow \mathbb{R}^{n\times n}$ where $g(u,v)$ is positive semi-definite for all $u,v\in\mathbb{R}$.

• Nice question. Did you try to go further when $F$ is uniform on $[0,1]$?
– Did
May 8, 2018 at 15:42
• Is there a reason why you don't want to use the Monte Carlo approach? May 8, 2018 at 15:48
• @Did I don't think the uniform case is actually any easier. The $F$ is not really the source of difficulty. I will have a play though.
– cfp
May 8, 2018 at 16:02
• @MikeHawk The Monte Carlo approach requires you to sample from a space of dimension $n$. While with $T$ samples the error is on the order of $\frac{1}{\sqrt{T}}$, in practice with a large $n$ (e.g. 2000), the constant is huge, and obtaining e.g. 1e-8 precision is essentially impossible.
– cfp
May 8, 2018 at 16:05
• Well, at least when $F$ is uniform, one knows every $E(X_{(k)})$, and, even though these do not appear in your computations, they seem to be needed to compute the covariance matrix.
– Did
May 8, 2018 at 16:22