I've been using simulation to investigate the behavior of cumulative density function evaluated at extreme values of a sample.

Functionally, my study has $n$ instrumented fish migrating downstream past a tracking station, according to an unknown probability density function with respect to time. With a sample size of $n$ fish, is it possible to calculate the expected value of the total probability density ($x\times100\%$ of the underlying distribution) represented by the min and max migration times in my sample? The ultimate goal is determining a sample size that will identify the time duration of $x\times100\%$ of the total downstream migration.

By means of brute-force simulation with n=10 and n=40 (and trying a few density functions)...

nsim <- 100000
x10 <- x40 <- list(NA, NA, NA, NA)

for(i in 1:nsim) {
  # normal, mean=0, sd=1
  x10[[1]][i] <- diff(pnorm(range(rnorm(10))))
  x40[[1]][i] <- diff(pnorm(range(rnorm(40))))

  # lognormal, meanlog=0, sdlog=1
  x10[[2]][i] <- diff(plnorm(range(rlnorm(10))))
  x40[[2]][i] <- diff(plnorm(range(rlnorm(40))))

  # exponential, rate=0
  x10[[3]][i] <- diff(pexp(range(rexp(10))))
  x40[[3]][i] <- diff(pexp(range(rexp(40))))

  # uniform, min=0, max=1
  x10[[4]][i] <- diff(punif(range(runif(10))))
  x40[[4]][i] <- diff(punif(range(runif(40))))

sapply(x10, mean)
## [1] 0.8181371 0.8182755 0.8181157 0.8182054

sapply(x40, mean)
## [1] 0.9512289 0.9512698 0.9511661 0.9512157

What I find truly interesting is that the total theoretical density doesn't seem to depend on the underlying density function! Not only that, but it seems to converge to $1-\frac{2}{n}$ (though not quite... some trial & error makes it look like $1-\frac{2}{n+1}$??)

Am I neglecting to think of something, or are there some interesting asymptotics at work?

Apologies for incorrect notation/terminology on my part - I've been out of grad school too long! Let me know if I can clarify anything.


"What I find truly interesting is that the total theoretical density doesn't seem to depend on the underlying density function!"

For continuous distributions this is not at all surprising. Consider transforming any continuous variate by its cdf $U=F_X(X)$ (i.e. the probability integral transform), which yields a (standard) uniform. When you do this, you don't change which point is the maximum, nor the height of the cdf there. So you can just work directly with the uniform, since any other continuous distribution will have the same property.

The next thing to notice is that $F_U(u)=u$, so for a standard uniform the height of the cdf is the argument itself. So now, you're just working out the expectation of the max and min. And for that, see Wikipedia on Order statistics of a uniform.

Now in your code, you appear to be looking at the range. $E(X_{(n)}-X_{(1)})$ = $E(X_{(n)})-E(X_{(1)})$ (properties of expectation), so you can just take the difference of the expectations. Now by symmetry $E(X_{(n)})=1-E(X_{(1)})$, so you just need $1-2E(X_{(1)})$. Now $E(X_{(1)})$ is the expectation of a $\text{Beta}(1,n)$ variate, which is $\frac{1}{n+1}$. Consequently the expected value of this "cdf-range" calculation for any continuous variate is $1-\frac{2}{n+1}$.

  • $\begingroup$ Thank you, Glen - that's elegant, clear, and helpful. $\endgroup$ – Matt Tyers Jun 28 '19 at 17:01
  • 1
    $\begingroup$ Note that these are exact small sample calculations; no asymptotic arguments were used. It works as well for n=2 as it does for n=10,000. $\endgroup$ – Glen_b Jun 29 '19 at 6:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.