# Sampling Order Statistics for Numerical Integration

This may be a stupid question. I want to do Monte Carlo integration over a region

$${\int}_{D_{1} \geq D_{2} \geq ... \geq D_{m} \geq 0} g(d_1,\ldots,d_m) f(d_1) f(d_2) \cdots f(d_m) \text{d}d_1\cdots\text{d}d_m$$

where $g$ is an arbitrary function and $f$ are density functions (either uniform or normal). My question is, for Monte Carlo integration over this region, can I simply draw $m$ samples from $f$ and sort them? If not, how do I sample points for this integration? If I were doing this in R, what function would I use to sample from this?

Thanks!

• Yes this is correct, with a m! correction to add [exercise] to compensate for the density of the ordered statistics. Simulation is thus straightforward if simulating from $f$ is straightforward. – Xi'an Nov 16 '14 at 21:21
• @Xi'an That looks pretty close to a decent answer. Please consider making it one -- I don't think you'd need to add much, since you cover the main issues in those two sentences. – Glen_b -Reinstate Monica Nov 16 '14 at 21:43
• steve if rf generates from density $f$, and g computes $g$, then something along the lines of replicate(nsim, {<...>} ) where <...> is a placeholder for a call to g with arguments being the sorted values generated from rf, the whole corrected by the factor for the ordering. – Glen_b -Reinstate Monica Nov 16 '14 at 21:53
• If this is for some subject, please add the self-study tag and read its tag wiki – Glen_b -Reinstate Monica Nov 16 '14 at 21:57
• Why do you refer to "$f$" in the plural? (If you intended to write $f_i(d_i)$, then the answer would be different.) – whuber Nov 17 '14 at 9:51

This is a correct assessment of the problem since $$\int_{D_{1} \geq D_{2} \geq ... \geq D_{m} \geq 0} g(d_1,\ldots,d_m) f(d_1) f(d_2) \cdots f(d_m) \text{d}d_1\cdots\text{d}d_m$$ can be turned into an expectation under the order statistics distribution $$m!\ \int g(d_1,\ldots,d_m) \dfrac{1}{m!}\mathbb{I}_{d_{1} \geq d_{2} \geq ... \geq d_{m} \geq 0} f(d_1) f(d_2) \cdots f(d_m) \text{d}d_1\cdots\text{d}d_m$$ (hence the missing factor $m!$). And simulating a sample from the order statistics is straightforward:
mean(g(t(apply(matrix(rt(1000,2),nr=100,nc=10),1,sort)))

where g is your function of interest (and rt is the $t$-distribution generator, to be replaced with yours).