Extreme value simulation with Monte Carlo

Question

I would like to seek your help with some questions to simulating extreme values. For example, I have written the following R code to generate QQplots for a normally distributed data, varying the size of the sample from 10 up to 1B. As seen in the graph as below, I had to go up to 1B samples to see extreme values @ 6 Sigma from the mean.

par(mfrow=c(2,3))
for(i in c(10, 100, 1e+3, 1e+4, 1e+5, 1e+6, 1e+7, 1e+8, 1e+9)){
  data <- rnorm(i, mean = 0, sd = 1)
  qqnorm(data, main=sprintf("Sample Size=%d", i)); qqline(data, col='red')
}

Question1: I am interested in extremes values only (tails). Though MC sampling generates most of the data close to the mean and I am wondering which method I can use to crate a few (up to a thousand maybe) samples and yet get some extreme values up to 6 sigma from mean.

Question2: the QQplot seems to be significant only after ~ 1K points. Also, no matter the sample size, the QQplot always shows some tail-offs. this must be some sampling error but I'm not sure. Would you please point me to any literature that explains this behavior ?

Answer 1

Suppose that $X$ follows a normal distribution with mean $\mu$ and variance $ \sigma^2$, that is $X \sim N(\mu, \sigma^2)$. The probability of getting a value at more that $6 \sigma$ from the mean is $$ \Pr\left(\left|\frac{X - \mu}{\sigma}\right|> 6\right) = \Phi(-6) + 1 - \Phi(6) = 2\Phi(-6) \approx 1.973175e-09 , $$ where the function $\Phi$ denotes the standard normal distribution function. That is about $1$ chance in $506'797'346$. So, no wonder that you need very large samples to observe such an event.

If you want a random sample from a normal distribution that contains such extreme events, then you'll have to generate very large samples, although you might occasionally observe them in smaller samples.

Now, if you're just interested in the tail(s) of the distribution, say at more than $4 \sigma$ from the mean, then you could sample from the tail of the distribution. To get a sample from the upper tail you could tailor the inverse transform sampling method to your needs, and sample uniform values in $[\Phi(4), 1]$, so $U \sim \mbox{unif}(\Phi(4), 1)$ and use $T = \mu + \sigma \Phi^{-1}(U)$. By construction, the random variable $T$ is at more than $4\sigma$ from the mean (in the upper tail). In fact, $T$ is a random variable from $N(\mu, \sigma^2)$ conditionned on the fact that $X > \mu + 4\sigma^2$.

Alternatively, you could rely on extreme value theory and use the Generalized Pareto distribution (GPD) to approximate of the tail of the normal distribution, and sample from the GPD.

As per the shape of the QQ-plots, one should expect (if the data come from the considered distribution) the points to be well alligned in larger samples but much less well alligned in smaller samples. In fact, this is closely linked to the asymptotic properties of the empirical distribution giving the pointwise convergence of the empirical distribution function to the true distribution function, which implies that empirical quantiles converge to true quantiles as the sample grows.