distribution for scaled Maximum of n independent Weibulls for $n \to \infty$

Question

Assume that $X_1, X_2,...\sim Weibull(\lambda, k) \quad iid.$, i.e.

$F(X_1\leq x) = 1-e^{-(\lambda x)^k}$

define $M_n:= \max\{X_1, ..., X_n\}$ and $\tilde{M}_n:=\frac{M_n-b_n}{a_n}$

according to extreme value theory the distribution of $\tilde{M}_n$ approaches for a proper choice of $a_n,b_n$ a distribution function:

$$G(x)= \exp\left\{-\left[1+\xi\left(\frac{x-\mu}{\sigma}\right)\right]^{-1/\xi}\right\}$$.

where $\mu \in \mathbb{R}, \sigma>0, \xi \in \mathbb{R}$ the case $\xi=0$ is interpreted as $\xi \to 0$

I.e.: $$\left\{1-\exp\left[-\lambda^k \left(\frac{x-b_n}{a_n}\right)^k\right]\right\}^n\to G(z), \quad n \to \infty$$

the expression on the left side can be reformulated as :

$$\left\{1-\exp\left[-\lambda^k \left(\frac{x-b_n}{a_n}\right)^k\right]\right\}^n=\left\{1- \frac{n\exp\left[-\lambda^k \left(\frac{x-b_n}{a_n}\right)^k\right]}{n}\right\}^n$$

so I need to find sequences $a_n,b_n$, such that

$$n\exp\left[-\lambda^k \left(\frac{x-b_n}{a_n}\right)^k\right]\to \left[1+\xi\left(\frac{x-\mu}{\sigma}\right)\right]^{-1/\xi}, n \to \infty$$

I would be glad if you could help me finding these sequences.

Answer 1

The determination of the domain of attraction and of the related constants $a_n$ and $b_n$ uses several functions related to the survival function, here given by $S(x) = \exp\{-(\lambda x)^k\}$ for $x > 0$. Of major importance are the tail-quantile function $U(t)$ and the hazard rate function $h(x)$.

The tail-quantile is obtained by solving $S(U) = 1/t$ for $U$, which leads to $U(t) = \lambda^{-1} (\log t)^{1/k}$ for $t> 1$. The hazard rate is $h(x) = -\text{d}\log S(x)/\text{d}x = \lambda^k \, k \,x^{k-1}$ for $x>0$. The derivative of the inverse hazard rate $1/h(x)$ is proportional to $x^{-k}$, hence tends to zero for large $x$ since $k>0$. So the Von Mises' condition holds. We know that the distribution is in the Gumbel domain of attraction and that we have the following convergence in distribution to the standard Gumbel

\begin{equation} \frac{M_n - U(n)}{a_n} \to \text{Gumbel}. \end{equation}

Moreover we know that we can choose $a_n$ as $e(U(n))$ where $e(x)$ denotes the mean residual life function given by $e(x) := A(x) / S(x)$ with $A(x) := \int_{x}^\infty S(t) \,\text{d}t$ for $x >0$.

We now proceed to the evaluation of $e(U(n))$ or to the determination of a quantity which is equivalent to it for large $n$. Using the change of variable $u := (\lambda t)^k$ we get

$$ A(x) = \int_{x}^\infty \exp\{-(\lambda t)^k\}\,\text{d}t= \frac{1}{\lambda k} \int_{(\lambda x)^k}^{\infty} u^{1/k -1} e^{-u} \text{d}u = (\lambda k)^{-1} \, \Gamma(s,\, v), $$

where $\Gamma(s,\,v)$ stands for the incomplete gamma function evaluated at $s:= 1/k$ and $v := (\lambda x)^k$. We can use the following known result about the incomplete gamma $\Gamma(s,\, v) \sim v^{s-1} e^{-v}$ for $v \to \infty$ which can be shown using integration by parts. So

$$ e(x) \sim \frac{(\lambda k)^{-1}\, \left\{ (\lambda x)^k \right\}^{1/k - 1} \exp\{- (\lambda x)^k \}}{S(x)} = \frac{x}{k} \, (\lambda x)^{-k}. $$

Note that $h(x) \times e(x)$ tends to $1$ for large $x$, which is clear from the last equivalence; this limit condition is both necessary an sufficient for the attraction to the Gumbel when $h(x)$ is monotonic for large $x$, as is the case here - see Theorem 1 in Galambos and Obretenov. We can choose $a_n$ as $1 / h(U(n))$, and our constants can be

$$ a_n = \dfrac{1}{\lambda k} \, (\log n)^{1/k -1}, \qquad b_n = \dfrac{1}{\lambda} \, (\log n )^{1/k}. $$

A precise statement of Von Mises' condition is found in the classical book Modelling Extremal Events by Embrechts P., Klüppelberg C. and Mikosch T. In this book (up to change in notations), the couple of constants is given in Table 3.4.4.

## Weibull parameters
k <- 2.5; lambda <- 10
## simulate
set.seed(123)
n <- 40; nsim <- 10000
X <- array(rweibull(n * nsim, shape = k, scale = 1/ lambda), dim = c(nsim, n))
M <- apply(X, 1, max)
bn <- log(n)^(1 / k) / lambda
an <- log(n)^(1 / k - 1) / lambda / k
Mscale <- scale(M, center = bn, scale = an)
hist(Mscale, breaks = 100, probability = TRUE, col = "lightyellow",
     main = sprintf("Simulated maxima for n = %d", n), xlab = "")
require(evd)
curve(dgumbel, add = TRUE, col = "orangered", lwd  = 2)