Numerické metody v R

$$ \xdef\scal#1#2{\langle #1, #2 \rangle} \xdef\norm#1{\left\lVert #1 \right\rVert} \xdef\dist{\rho} \xdef\and{\&}\xdef\brackets#1{\left\{ #1 \right\}} \xdef\parc#1#2{\frac {\partial #1}{\partial #2}} \xdef\mtr#1{\begin{pmatrix}#1\end{pmatrix}} \xdef\bm#1{\boldsymbol{#1}} \xdef\mcal#1{\mathcal{#1}} \xdef\vv#1{\mathbf{#1}}\xdef\vvp#1{\pmb{#1}} \xdef\ve{\varepsilon} \xdef\l{\lambda} \xdef\th{\vartheta} \xdef\a{\alpha} \xdef\tagged#1{(\text{#1})} \xdef\tagged*#1{\text{#1}} \xdef\tagEqHere#1#2{\href{#2\#eq-#1}{(\text{#1})}} \xdef\tagDeHere#1#2{\href{#2\#de-#1}{\text{#1}}} \xdef\tagEq#1{\href{\#eq-#1}{(\text{#1})}} \xdef\tagDe#1{\href{\#de-#1}{\text{#1}}} \xdef\T#1{\htmlId{eq-#1}{#1}} \xdef\D#1{\htmlId{de-#1}{\vv{#1}}} \xdef\conv#1{\mathrm{conv}\, #1} \xdef\cone#1{\mathrm{cone}\, #1} \xdef\aff#1{\mathrm{aff}\, #1} \xdef\lin#1{\mathrm{Lin}\, #1} \xdef\span#1{\mathrm{span}\, #1} \xdef\O{\mathcal O} \xdef\ri#1{\mathrm{ri}\, #1} \xdef\rd#1{\mathrm{r}\partial\, #1} \xdef\interior#1{\mathrm{int}\, #1} \xdef\proj{\Pi} \xdef\epi#1{\mathrm{epi}\, #1} \xdef\grad#1{\mathrm{grad}\, #1} \xdef\hess#1{\nabla^2\, #1} \xdef\subdif#1{\partial #1} \xdef\co#1{\mathrm{co}\, #1} $$

Rychlost konvergence

Definice $\D{3.1}$

Nechť jsou dány 2 posloupnosti $\brackets{e_ k}_ {k = 0}^\infty$ a $\brackets{h_ k}_ {k = 0}^\infty$ takové, že $$ e_k \in [0, \infty), \quad e_k \to 0 \quad \and \quad h_k \in [0, \infty), \quad h_k \to 0. $$ Řekneme, že posloupnost $\brackets{e_k}$ konverguje rychleji (pomaleji) než $\brackets{h_k}$, pokud existuje index $\tilde k \in \N_0$ takový, že $$ e_k \leq_{(\geq)} h_k \quad \forall k \in [\tilde k, \infty) \cap \N_0 $$

Definice $\D{3.2}$ (Rychlost konvergence)

Nechť je dána posloupnost $\brackets{e_ k}_ {k = 0}^\infty$ splňující $e_k \in [0, \infty)$ a $e_k \to 0$. Řekneme, že posloupnost $\brackets{e_k}$ konverguje

alespoň lineárně s rychlostí $\beta \in (0,1)$, pokud konverguje rychleji než geometrická posloupnost se členy tvaru $q \bar \beta^k$, kde $q > 0$ a $\bar \beta \in (\beta, 1)$.
- čím větší $\bar \beta$, tím pomaleji jde tato geometrická posloupnost k nule - tj. konverguje rychleji než geometrická posloupnost s $\bar \beta$ větší než $\beta$
nejvýše lineárně s rychlostí $\beta \in (0,1)$, pokud konverguje pomaleji než geometrická posloupnost se členy tvaru $q \bar \beta^k$, kde $q > 0$ a $\bar \beta \in (0, \beta)$.
lineárně s rychlostí $\beta \in (0,1)$, pokud konverguje nejvýše a současně alespoň lineárně s rychlostí $\beta$.
superlineárně (sublineárně), pokud konverguje rycheji (pomaleji) než libovolná geometrická posloupnost se členy tvaru $q \beta^k$, kde $q > 0$ a $\beta \in (0,1)$.

Definice $\D{3.4}$

Nechť je dána posloupnost $\brackets{e_ k}_ {k = 0}^\infty$ splňující $e_k \in [0, \infty)$ a $e_k \to 0$, přičemž $\brackets{e_k}$ konverguje superlineárně. Řekneme, že posloupnost $\brackets{e_k}$ konverguje

alespoň superlineárně s řádem $p > 1$, pokud konverguje rychleji než všechny posloupnosti se členy tvaru $q \beta^{\bar p^k}$, kde $q > 0$ a $\beta \in (0,1)$ a $\bar p \in (1, p)$
- čím větší $\bar p$, tím rychleji posloupnost $q \beta^{\bar p^k}$ konverguje - tj. $\brackets{e_k}$ konverguje rychleji než všechny posloupnosti s menším $p$
nejvýše superlineárně s řádem $p > 1$, pokud konverguje pomaleji než všechny posloupnosti se členy tvaru $q \beta^{\bar p^k}$, kde $q > 0$ a $\beta \in (0,1)$ a $\bar p \in (p, \infty)$
superlineárně s řádem $p > 1$, pokud konverguje nejvýše a současně alespoň superlineárně s řádem $p$.
superlineárně s řádem $p = 1$, pokud konverguje pomaleji než všechny posloupnosti se členy tvaru $q \beta^{\bar p^k}$, kde $q > 0$ a $\beta \in (0,1)$ a $\bar p \in (1, \infty)$

Metody

Definice $\D{3.1.1}$ (Unimodální funkce)

Nechť je dán interval $I \subset \R$ a funkce $f: I \to \R$. Řekneme, že $f$ je unimodální na $I$, jestliže existuje $x^* \in I$ takové, že

$f(x_1) > f(x_2)$ pro libovolná $x_1, x_2 \in I$ splňující $x^* > x_1 > x_2$
$f(x_1) < f(x_2)$ pro libovolná $x_1, x_2 \in I$ splňující $x^* < x_1 < x_2$

Jinými slovy, unimodální funkce je klesající na $(-\infty, x^* ) \cap I$ (tj. nalevo od $x^* $) a rostoucí na $I \cap (x^* , \infty)$ (tj. napravo od $x^* $).

Unimodalita neimplikuje konvexnost (ani se spojitostí), pouze kvazikonvexnost

Naopak, konvexní funkce nemusí nutně být unimodální (ale ostrá konvexnost $\implies$ unimodalita)

Konvexní funkce nemusí být např. jen rostoucí, ale i neklesající

V této části budeme řešit úlohu $$ f(x) \to \min, \qquad x \in I := [a,b] \tag{\T{3.1.1}} $$

Lemma $\D{3.1.2}$

Nechť $f: I \to \R$ je unimodální na $I$ a $x_1, x_2 \in I$ jsou takové, že $x_1 < x_2$.

Je-li $f(x_1) \leq f(x_2)$, pak $x^* \leq x_2$
Je-li $f(x_1) \geq f(x_2)$, pak $x^* \geq x_1$

Upozornění:
Dále uvažujme POUZE UNIMODÁLNÍ funkce.

Navíc nechť $N$ značí povolený počet vyčíslení a přesnost těchto metod je dáno jako $|\bar x - x^* |$, kde $x^* $ je přesné řešení úlohy $\tagEq{3.1.1}$ a $\bar x$ jeho nalezená aproximace.

Metoda prostého dělení

Tato metoda není efektivní a je to de facto hrubá síla

Podle parity $N$ určíme dělící body intervalu $I$.

$N$ liché	$N$ sudé
$$x_i := a + {b - a \over N + 1} i, \quad i=1, \dots, N = 2k - 1$$	$$x_{2i} := a + {b - a \over k + 1} i \quad \and \quad x_{2i - 1} := x_{2i} - \delta, \ i = 1, \dots, k := N/2,$$

kde $\delta$ je vhodné malé číslo.

Poté vyčíslíme $f(x_1), \dots, f(x_N)$ (což v případě $N$ sudého a $\delta \in \set{0, {b-a \over k+ 1}}$ znamená pouze $k$ vyčíslení) a nechť v $x_j$ nastává nejmenší hodnota, tj. $$ f(x_j) = \min_{1 \leq i \leq N} f(x_i) $$ Pak z Lemma $\tagDe{3.1.2}$ plyne, že $x^* \in [x_{j-1}, x_{j+1}]$ a tento interval nazveme interval lokalizace minima (ILM) a za aproximaci $x^* $ vezmeme střed ILM, tj. $\bar x := {x_{j-1} + x_{j+1} \over 2}$

Pro délku $l_N$ intervalu lokalizace minima platí $$ l_N := \max_{1 \leq i \leq N} (x_{i+1} - x_{i-1}) = \begin{cases} 2 {b - a \over N+1}, & N = 2k - 1 \ {b - a \over (N/2) + 1} + \delta, & N = 2k \end{cases} $$

Pro $N$ sudé je poslední interval delší, proto dostáváme takový tvar $l_N$

Přesnost této metody je dána polovinou ILM, tj. ${l_N \over 2}$

Rychlost konvergence této metody je sublineární, navíc je tento algoritmus pasivní, tj. volba $x_{m+1}$ nezáleží na $x_1, \dots, x_m$ (závisí pouze na $N$, či na $N$ a $\delta$).

Metoda půlení intervalu

Nechť nyní $N = 2k$. Položme $a_0 = a$, $b_0 = b$ a $$ x_i^- := {a_{i-1} + b_{i-1} \over 2} - \delta \quad \and \quad x_i^+ := {a_{i-1}+b_{i-1} \over 2} + \delta, $$ kde $\delta > 0$ je dostatečně malé a $i = 1, \dots, k$. Vyčíslíme funkci v $x_i^-, x_i^+$, tj. dostaneme $f(x_i^-), f(x_i^+)$. Pak

jestliže $f(x_i^-) < f(x_i^+)$, pak podle Lemma $\tagDe{3.1.2}$ je ILM $[a_{i-1}, x_i^+] \implies a_i = a_{i-1}, b_i = x_i^+$
jestliže $f(x_i^-) > f(x_i^+)$, pak podle Lemma $\tagDe{3.1.2}$ je ILM $[x_i^-, b_{i-1}] \implies a_i = x_i^-, b_i = b_{i-1}$

Takto můžeme tento proces opakovat ($k$-krát, jelikož máme $N = 2k$ povolených vyčíslení), kdy za $a,b$ volíme krajní body ILM pro každý krok. Zřejmě, jako aproximaci $x^* $ v $k$-tém kroku bereme střed ILM pro $k$-tý krok.

Délka ILM je v tomto případě $$ l_k = {b - a\over 2^k} + {(2^k - 1)\delta \over 2^{k-1}}, $$ přičemž $\delta \in [0, {b-a\over 2}]$ a navíc pro $k \to \infty$ je $l_k \to 2\delta$.

Z tohoto vyplývá, že čím menší $\delta$, tím je metoda přesnější. Nicméně ve skutečnosti se můžeme dostat k zaokrouhlovacím chybám, které dokonce mohou způsobit, že špatně určíme velikosti $f(x_i^-), f(x_i^+)$ (tím pádem bychom řekli, že $x^* $ je v opačném intervalu, než je ve skutečnosti)

Tato metoda konverguje lineárně s rychlostí $1/2$.

Metoda zlatého řezu

Myšlenka metody zlatého řezu "vylepšuje" metodu půlení intervalu tak, že každá další iterace umožňuje pouze jedno další vyčíslení.

Zde $\tau$ je řešení rovnice $\tau^2 - \tau - 1 = 0$, tj. a $\frac 1 \tau \approx 0.618$

Mějme funkci $f$, interval $[a,b]$, přesnost $\ve$ nebo počet vyčíslení $N \geq 2$:

(Inicializace) Položíme $a_0 := a, b_0 := b$ a $k := 1$. Vypočteme
$$\l_1 := a_0 + {b_0 - a_0 \over \tau^2} \quad \and \quad \mu_1 := a_0 + {b_0 - a_0 \over \tau}$$
Je-li $k = N$, pokračujeme částí 5., jinak následuje krok 3.
Vyčíslíme $f(\l_k)$ a $f(\mu_k)$. Jestliže $f(\l_k) \geq f(\mu_k)$:
1. Položíme $a_k := \l_k, b_k = b_{k-1}, \l_{k+1} := \mu_k$ a
  $$f(\l_{k+1}) := f(\mu_k), \quad \mu_{k+1} := a_k + {b_k - a_k \over \tau}$$ a pokračujeme na krok 4.
2. Položíme $a_k := a_{k-1}, b_k := \mu_k, \mu_{k+1} := \l_k$ a
  $$f(\mu_{k+1}) := f(\l_k), \quad \l_{k+1} := a_k + {b_k - a_k \over \tau^2}$$ a pokračujeme na krok 4.
Položíme $k := k+1$ a pokračujeme krokem 2.
Stanovíme poslední ILM jako $[a_{k-1}, b_{k-1}]$ a vypočteme $\bar x := {a_{k-1} + b_{k-1} \over 2}$. KONEC

Tato metoda konverguje lineárně s rychlostí $\frac 1 \tau \approx 0.618$

Toto neznamená, že by na stejný počet vyčíslení byla tato metoda horší než metoda půlení intervalu

Fibonacciho metoda

V této poslední metodě uvažujme, že zkrácení $\delta$ může být jiné v každé kroku metody.

Nechť $F_n$ je $n$-té Fibonacciho číslo

Máme povoleno $N$ vyčíslení, takže $M = N - 1$ a $$ \l_i = a_{i - 1} + {F_{N - i -1} \over F_{N - i + 1}} l_{i-1} = b_{i-1} - {F_{N - 1} \over F_{N - 1 + i}} l_{i-1} $$ $$ \mu_i = a_{i-1} + {F_{N - 1} \over F_{N - 1 + i}} l_{i-1} = b_{i-1} - {F_{N - i -1} \over F_{N - i + 1}} l_{i-1} $$

Tato metoda konverguje lineárně s rychlostí $\frac 1 \tau \approx 0.618$, tj. stejně jako metoda zlatého řezu