12  Rozdělení extrémních hodnot

V hydrologii a klimatologii se často nezajímáme o „běžné“ hodnoty náhodné proměnné, ale o její extrémy – například maximální roční průtoky, nejvyšší denní úhrny srážek nebo minimální průtoky v období sucha. K jejich popisu se používá teorie rozdělení extrémních hodnot (Extreme Value Theory, EVT).

Základní myšlenkou je studium chování maxim (případně minim) náhodného výběru místo samotné distribuční funkce původních dat.

12.1 Vymezení extrémní hodnoty

12.1.1 Bloková maxima

12.1.2 Překročení limitu

Z dat jsou užita pouze pozorování překračující stanovenou mez, například 95 % kvantil.

Úloha

1. 1. Nahrajte do prostředí data z klementiské řady měření
   2. Metodou blokových maxim stanovte maxima

klementinum_short <- klementinum |> 
  dplyr::filter(yr %in% 2000:2020) 
window <- 2000:2020


i <- seq(from = as.Date("2000-01-01"), to = as.Date("2020-12-31"), by = "1 year")
j <- ifelse(as.numeric(i) %% 2 == 0, "black", "#00000090")

plot(x = klementinum_short$date,
     y = klementinum_short$tavg, 
     type = "p", 
     cex = 0.1, 
     col = "#00000050", 
     ylab = expression(italic(T)[avg]*degree*C), 
     xlab = "")
## draw rectangles with bottom left (100, 300)+i
## and top right (150, 380)+i
rect(0+i, 365+i, -10, 20, col = j, 
     border = j)
abline(h = threshold, 
       col = "red", 
       lwd = 2)

Klementinum 2000-2020. Identifikace extrémů s pomocí ‘Peak-over-threshold’.

12.2 Gumbellovo rozdělení

Hustota funkce Gumbellova rozdělení je dána předpisem

\[ f(x) = \dfrac{1}{d}\exp\left(-\dfrac{x-c}{d}\right)\cdot\exp\left[-\exp\left(-\dfrac{x-c}{d}\right)\right] \] kde \(c\) a \(d\) jsou parametry. Kvantilová funkce je pak \[ x_T = c - d\cdot\ln\left[-\ln\left(1-\dfrac{1}{T}\right)\right] \]

kde \(T\) je doba opakování a parametry odhadnuté metodou mometů jsou:

\[ d = \dfrac{\sqrt{6}}{\pi}\sigma, \quad c = \mu-0.5772\cdot d \]

a \(\mu\) \(\sigma\) jsou momenty celkového souboru.

fitgumbel <- function(x, ...) {
  mX <- mean(x)
  sdX <- sd(x)
  dX <- sqrt(6)/pi*sdX
  cX <- mX - 0.5772*dX
  curve(cX - d*log(-log(1-1/x)), add = TRUE, ...)
}

plot(x = c(1,100), 
     y = c(0,100), 
     type = "n", 
     xlab = "Doba opakování (roky)", 
     ylab = "Maximální denní průtok v roce (mm/d)")

# fitgumbel()