Kriging
Author
Albert FloresKriging či počeštěně krigování je ve statistice (původně v geostatistice) metoda interpolace, kde jsou interpolované hodnoty modelovány gaussovským procesem podle apriorních kovariancí. Za vhodných předpokladů dává kriging nejlepší lineární nestrannou předpověď střední hodnoty. Interpolační metody založené na jiných kritériích, jako je například hladkost, nemusejí přinést nejpravděpodobnější střední hodnoty. Tato metoda se běžně používá v oblasti prostorové analýzy a počítačových experimentů. Technika je také známá jako Kolmogorova-Wienerova predikce.
V roce 1960 teoretický základ metody vypracoval francouzský matematik Georges Matheron na základě diplomové práce Danieho G. +more Krigeho, průkopníka ve vykreslování vzdáleností vážených průměrných hodnot zlata v útesovém komplexu Witwatersrand v Jihoarfické republice. Krige se snažil odhadnout nejpravděpodobnější rozložení zlata na základě vzorků z několika vrtů.
Hlavní zásady
Související pojmy a techniky
Jednoduchá myšlenka je předpoklad funkční hodnoty daného bodu vypočítaného váženého průměru známých funkčních hodnot sousedních bodů. Metoda je matematicky úzce související s regresní analýzou. +more Obě teorie odvozují nejlepší lineární nestranný odhad, založený na předpokladu kovariance. Využitím gauss-Markovy teorie dokážeme nezávislost odhadů a omylů a využívá velmi podobné vzorce. Nicméně jsou použitelné rozdílné rámce: kriging je dělaný pro odhad jednu realizaci z náhodných polí, zatímco regresní model je založen na základě pozorování vícerozměrné datové sady.
Odhad krigingu může být také viděn jako křivka v reprodukujícím se jádra Hilbertova prostoru s reprodukujícím jádrem dané kovarianční funkce. Na rozdíl od klasického krigingu je provedena interpretace: zatímco křivka je motivována minimální normou interpolací založenou na Hilbertově prostorové struktuře, kriging je motivován očekávanou kvadratickou chybou na základě stochastického modelu.
Kriging s polynomiálním trendem ploch je matematicky identický jako zobecněný polygon nejmenších čtverců křivky.
Kriging může být známý jako forma Bayesianova závěru. Kriging začíná s základním rozšířením přes funkce. +more Tento základ bere formu gausova procesu. N vzorků z funkce která má normální rozdělení a kovariance mezi dvěma vzorky je kovariance funkce (nebo jádra) gausova procesu hodnoceny prostorovým rozmístěním dvou bodů. Sada hodnot je pak pozorována, každá hodnota souvisí s prostorovým rozmístěním. Takže nová hodnota může být předpovězena v nějakém novém prostorovém rozmístění, kombinací základní gausovskou funkcí s gausovskou věrohodnostní funkcí pro každou pozorovanou hodnotu. Výsledná distribuce pozadí je také Gaussian s průměrnou hodnotou a kovariancí, která může být snadno vypočítána z pozorovaných hodnot, jejich rozdílnosti a jádra matic odvozených od předchozích.
Geostatistický odhad
V geostatistickém modelu jsou některá data interpretována jako výsledky náhodného procesu. Skutečnost, že tyto modely obsahují nejistotu v jeho konceptualizaci, neznamená, že tento jev - les, zavodnění nebo minerální usazenina - mají výsledky z náhodného procesu ale pouze umožňují vytvořit metodologické základy pro prostorové závěry z množství nevypozorovaných lokalit a kvantifikovat nejistoty spojené s odhadem. +more Stochastický proces je v rámci tohoto modelu, způsob jak přistupovat k souboru údajů shromážděných ze vzorků. První krok v geostatistických modulacích je vytvoření náhodného procesu, který nejlépe popisuje sadu experimentálních pozorovaných dat. Hodnota prostorově umístěná v x_1 (obecné označení datové sady s geografickým souřadnicovým systémem) je interpretovaná jako realizace z(x_1) náhodné proměnné Z(x_1). V prostoru A, kde datová sada vzorků je rozptýlená, existuje N realizace těchto náhodných veličinZ(x_1), Z(x_2), \cdots, Z(x_N), korelujících mezi sebou. Datová sada náhodných veličin, představuje náhodnou funkci pouze u těch u kterých je známá pouze jedna realizace z(x_i) - soubor experimentálních dat. S pouze jednou realizací každé náhodné proměnné, je teoreticky možné určit všechny statistické parametry jednotlivých proměnných nebo funkcí. :Navrhované řešení geostatistických formalismů spočívá v předpokladu různých stupňů stacionarity v náhodné funkci, aby bylo možné se domnívat závěru některých statistických hodnot.
Například, jestliže skupina vědců předpokládá odpovídající, je vhodné na základě homogenity vzorků v oblasti A, kde je distribuována proměnná, hypotetizovat, že první část je stacionární (všechny náhodné proměnné mají stejný průměr), z čehož vyplývá, že průměr může být odhadnut aritmetickým průměrem vybraných hodnot. Soudě hypotézu jako je tato, tak v případě potřeby dostatečně ověřte homogenitu vzorových hodnot, jestli jsou reprezentativní.
Hypotéza stacionarity týkající se druhé části je definována následujícím způsobem: korelace mezi dvěma náhodnými proměnnými závisí pouze na prostorové vzdálenosti, která je odděluje a je závislá na jeho poloze:
C(Z(x_1),Z(x_2)) = C(Z(x_i),Z(x_i+\mathbf{h})) = C(\mathbf{h})
\gamma(Z(x_1),Z(x_2)) = \gamma(Z(x_i),Z(x_i+\mathbf{h})) = \gamma(\mathbf{h})
kde \mathbf{h} = (x_1,x_2) = (x_i,x_i+\mathbf{h})
Tato hypotéza umožňuje odvozovat kroky - variogram a kovariogram - na základě N vzorků:
\gamma(\mathbf{h})=\frac{1}{2N(\mathbf{h})}\sum^{N(\mathbf{h})}_{i=1}\left(Z(x_i)-Z(x_i+\mathbf{h})\right)^2
C(\mathbf{h})=\frac{1}{N(\mathbf{h})}\sum^{N(\mathbf{h})}_{i=1}\left(Z(x_i)Z(x_i+\mathbf{h})\right)-m(x_i)m(x_i+\mathbf{h})
kde m(x_i)=\frac{1}{N(\mathbf{h})}\sum^{N(\mathbf{h})}_{i=1}Z(x_i)
Lineární odhad
Prostorový závěr nebo odhad kvantityZ: \mathbb{R}^n\rightarrow\mathbb{R} na neurčeném místě x_0 je vypočítán z lineární kombinace pozorovaných hodnotz_i=Z(x_i) a vah w_i(x_0),\;i=1,\ldots,N:
\hat{Z}(x_0) = \begin{bmatrix} w_1 & w_2 & \cdots & w_N \end{bmatrix} \cdot \begin{bmatrix} z_1\\ z_2\\ \vdots\\ z_N \end{bmatrix} = \sum_{i=1}^n w_i(x_0) \times Z(x_i)
Váhy w_ijsou určeny ke shrnutí dvou extrémně důležitých procedur v závěru prostorového procesu: * odrážet strukturální "blízkost" vzorků odhadovaného umístění, x_0 * ve stejnou dobu by měli mít degradující efekt, aby se předešlo zkreslení případnými shluky
Když vypočítáváme váhy w_i, tak jsou v geostatistice dva závěry: nestranný a minimální odchylka odhadu.
Jestliže mračno reálných hodnot Z(x_0)je vykresleno proti odhadovaným hodnotám \hat{Z}(x_0)tak kritérium pro souhrnnou nestrannost, vnitřní stacionarita nebo široký smysl stacionarity oboru, znamená, že průměr odhadů musí být roven průměru skutečných hodnot.
Druhé kritérium říká, že průměrná kvadratická odchylka(\hat{Z}(x)-Z(x)) musí být minimální, což znamená, že když se mračno odhadovaných hodnot oproti mraku skutečných hodnot více rozptýlí, odhad je více nepřesný.
Metody
V závislosti na stochastických vlastnostech náhodného pole a různým stupněm předpokládané stacionarity, mohou vzniknout různé metody pro výpočet hmotnosti. Klasické metody jsou:
* Běžný kriging předpokládá stacionaritu prvního momentu všech náhodných veličin: E\{Z(x_i)\}=E\{Z(x_0)\}=m, kde m je neznámé. * Jednoduchý kriging předpokládá známý stacionární průměr: E\{Z(x)\}=m, kde m je známé. +more * Univerzální kriging počítá s obecným polynomickým modelem trendu, např. lineární model trendu E\{Z(x)\}=\sum_{k=0}^p \beta_k f_k(x). * IRFk-kriging předpokládá, že E\{Z(x)\} je polynomická neznámá v x. * Indikační kriging používá indikační funkce namísto samotného procesu, aby bylo možné odhadnout pravděpodobnosti přechodů. * Více-indikační kriging je verze indikačního krigingu pracující s rodinou ukazatelů. Nicméně, V-IK přestal v posledních letech vyhovovat jako interpolační technika. K tomu došlo hlavně díky problémům spojeným s provozem a validací modelu. Podmíněná simulace se v tomto oboru rychle stává uznávanou náhradní technikou. * Disjunktivní kriging je nelineární zobecnění krigingu. * Logaritmicko-normální kriging interpoluje pozitivní data pomocí logaritmů.
Běžný kriging
Neznámá hodnota Z(x_0) je interpretována jako náhodná proměnná nacházející se v x_0, jakož i hodnoty sousedních vzorků Z(x_i), i=1,\cdots ,N. Odhad \hat{Z}(x_0) je také interpretován jako náhodná proměnná nacházející se v x_0, výsledek lineární kombinace proměnných.
Aby bylo možné vyvodit krigovací systém pro předpoklad modelu, následující chyba při odhadování Z(x) v x_0 je deklarována:
:\epsilon(x_0) = \hat{Z}(x_0) - Z(x_0) = \begin{bmatrix}W^T&-1\end{bmatrix} \cdot \begin{bmatrix}Z(x_i)&\cdots&Z(x_N)&Z(x_0)\end{bmatrix}^T = \sum^{N}_{i=1}w_i(x_0) \times Z(x_i) - Z(x_0)
Obě kritéria kvality uvedená dříve mohou nyní být vyjádřena z hlediska střední hodnoty a rozptylu nové náhodné proměnné \epsilon(x_0):
Objektivnost
Vzhledem k tomu, že náhodná funkce je stacionární E(Z(x_i))=E(Z(x_0))=m, uvádí se následující omezení:
:E\left(\epsilon(x_0)\right)=0 \Leftrightarrow \sum^{N}_{i=1}w_i(x_0) \times E(Z(x_i)) - E(Z(x_0))=0 \Leftrightarrow : :\Leftrightarrow m\sum^{N}_{i=1}w_i(x_0) - m=0 \Leftrightarrow \sum^{N}_{i=1}w_i(x_0) = 1 \Leftrightarrow \mathbf{1}^T \cdot W = 1
Aby se zajistilo, že model je objektivní, součet vah musí být jedna.
Minimální variance: minimalizuje E\left(\epsilon(x_0)^2\right)
Dva odhady mohou mít \epsilon(x_0)=0, ale disperze kolem jejich střední hodnoty určí rozdíl mezi kvalitou odhadů.
:\begin{array}{rl} Var(\epsilon(x_0)) &= Var\left(\begin{bmatrix}W^T&-1\end{bmatrix} \cdot \begin{bmatrix}Z(x_i)&\cdots&Z(x_N)&Z(x_0)\end{bmatrix}^T\right) =\\ &\overset{*}{=} \begin{bmatrix}W^T&-1\end{bmatrix} \cdot Var\left(\begin{bmatrix}Z(x_i)&\cdots&Z(x_N)&Z(x_0)\end{bmatrix}^T\right) \cdot \begin{bmatrix}W\\-1\end{bmatrix} \end{array}
* navštivte covariance matrix pro detailní vysvětlení (v angličtině)
:Var(\epsilon(x_0)) \overset{*}{=} \begin{bmatrix}W^T&-1\end{bmatrix} \cdot \begin{bmatrix}Var_{x_i}& Cov_{x_ix_0}\\Cov_{x_ix_0}^T & Var_{x_0}\end{bmatrix} \cdot \begin{bmatrix}W\\-1\end{bmatrix}
* kde literály představují \left\{Var_{x_i}, Var_{x_0}, Cov_{x_ix_0}\right\} stand for \left\{Var\left(\begin{bmatrix}Z(x_1)&\cdots&Z(x_N)\end{bmatrix}^T\right), Var\left(Z(x_0)\right), Cov\left(\begin{bmatrix}Z(x_1)&\cdots&Z(x_N)\end{bmatrix}^T,Z(x_0)\right)\right\}.
Jednou definovaná kovariance modelu nebo variogram, C(\mathbf{h}) nebo \gamma(\mathbf{h}), platné ve všech oblastech analýzy Z(x), pak můžeme definovat vztah pro odhad rozptylu libovolného odhadu v závislosti na kovarianci mezi vzorky:
:\left\{\begin{array}{l} Var(\epsilon(x_0)) = W^T \cdot Var_{x_i} \cdot W - Cov_{x_ix_0}^T \cdot W - W^T \cdot Cov_{x_ix_0} + Var_{x_0}\\ Var(\epsilon(x_0)) = Cov(0) + \sum_{i}\sum_{j}w_iw_jCov(x_i,x_j) - 2 \sum_iw_iC(x_i,x_0)\end{array} \right.
Některé závěry mohou být z těchto výrazů stanoveny. Rozptyl odhadu:
* není kvantifikovatelný na kterýkoliv lineární odhad, jakmile se předpokládá stacionarita průměru a prostorových kovariancí nebo variogramy. * roste, když se kovariance mezi vzorky a místem pro odhad snižuje. +more To znamená, že když vzorky jsou daleko od x_0, nejhorší je odhad. * roste zároveň s "prioritním" rozptylem C(0) proměnné Z(x). Pokud je proměnná méně rozptýlená, rozptyl je nižší v každém bodě prostoru A. * nezávisí na hodnotách vzorků. To znamená, že stejné prostorové uspořádání (se stejnými geometrickými vztahy mezi vzorky a body odhadu) vždy reprodukuje stejný odhad rozptylu v jakékoli části areálu A. Tímto způsobem rozptyl neměří nejistotu odhadu produkovanou místní proměnnou.
Systém Kriging
:\begin{align} &\underset{W}{\operatorname{minimize}}& & W^T \cdot Var_{x_i} \cdot W - Cov_{x_ix_0}^T \cdot W - W^T \cdot Cov_{x_ix_0} + Var_{x_0} \\ &\operatorname{subject\;to} & &\mathbf{1}^T \cdot W = 1 \end{align}
Řešení tohoto optimalizačního problému má za následek vznik systému kriging:
:\begin{bmatrix}\hat{W}\\\mu\end{bmatrix} = \begin{bmatrix} Var_{x_i}& \mathbf{1}\\ \mathbf{1}^T& 0 \end{bmatrix}^{-1}\cdot \begin{bmatrix}Cov_{x_ix_0}\\ 1\end{bmatrix} = \begin{bmatrix} \gamma(x_1,x_1) & \cdots & \gamma(x_1,x_n) &1 \\ \vdots & \ddots & \vdots & \vdots \\ \gamma(x_n,x_1) & \cdots & \gamma(x_n,x_n) & 1 \\ 1 &\cdots& 1 & 0 \end{bmatrix}^{-1} \begin{bmatrix}\gamma(x_1,x^*) \\ \vdots \\ \gamma(x_n,x^*) \\ 1\end{bmatrix}
další parametr \mu je Lagrangeův násobič použitý při minimalizaci chyby krigingu \sigma_k^2(x) pro dosažení nestrannosti.
Jednoduchý kriging
Jednoduchý kriging je matematicky nejjednodušší, ale nejméně obecný. To předpokládá, že očekávaná hodnota náhodného pole je známá, a opírá se o kovarianční funkci. +more Nicméně, ve většině aplikací nejsou očekávání, ani kovariance známy předem.
Praktické předpoklady pro uplatnění jednoduchého krigingu jsou: * Široká stacionarita pole. * Očekávání je všude nulové: \mu(x)=0. +more * Známá kovarianční funkce c(x,y)=\mathrm{Cov}(Z(x),Z(y)).
Systém Kriging
Krigovací váhy jednoduchého krigingu nemají nestrannou podmínku a jsou dány systémem rovnic jednoduchého krigingu: :\begin{pmatrix}w_1 \\ \vdots \\ w_n \end{pmatrix}= \begin{pmatrix}c(x_1,x_1) & \cdots & c(x_1,x_n) \\ \vdots & \ddots & \vdots \\ c(x_n,x_1) & \cdots & c(x_n,x_n) \end{pmatrix}^{-1} \begin{pmatrix}c(x_1,x_0) \\ \vdots \\ c(x_n,x_0) \end{pmatrix}
Toto je analogické s lineární regresí Z(x_0) na jiném z_1 , \ldots, z_n.
Odhad
Interpolace pomocí jednoduchého krigingu je dána: :\hat{Z}(x_0)=\begin{pmatrix}z_1 \\ \vdots \\ z_n \end{pmatrix}' \begin{pmatrix}c(x_1,x_1) & \cdots & c(x_1,x_n) \\ \vdots & \ddots & \vdots \\ c(x_n,x_1) & \cdots & c(x_n,x_n) \end{pmatrix}^{-1} \begin{pmatrix}c(x_1,x_0) \\ \vdots \\ c(x_n,x_0)\end{pmatrix}
Krigovací chyba ja dána: :\mathrm{Var}\left(\hat{Z}(x_0)-Z(x_0)\right)=\underbrace{c(x_0,x_0)}_{\mathrm{Var}(Z(x_0))}- \underbrace{\begin{pmatrix}c(x_1,x_0) \\ \vdots \\ c(x_n,x_0)\end{pmatrix}' \begin{pmatrix} c(x_1,x_1) & \cdots & c(x_1,x_n) \\ \vdots & \ddots & \vdots \\ c(x_n,x_1) & \cdots & c(x_n,x_n) \end{pmatrix}^{-1} \begin{pmatrix}c(x_1,x_0) \\ \vdots \\ c(x_n,x_0) \end{pmatrix}}_{\mathrm{Var}(\hat{Z}(x_0))}
což vede ke generalizované verzi nejmenších čtverců Gauss-Markovova teorému: :\mathrm{Var}(Z(x_0))=\mathrm{Var}(\hat{Z}(x_0))+\mathrm{Var}\left(\hat{Z}(x_0)-Z(x_0)\right).
Vlastnosti
Odhad krigingu je nestranný:: E[\hat{Z}(x_i)]=E[Z(x_i)] * Odhad krigingu využívá skutečně pozorované hodnoty: \hat{Z}(x_i)=Z(x_i) (za předpokladu, že nevznikla chyba ve měření, při kterém vznikly) * Odhad krigingu \hat{Z}(x) je nejlepší lineární nestranný odhad pro Z(x) v případě, že předpoklad trvá. Avšak: ** Stejně jako u jiných metod: Pokud předpoklady netrvají, krigování může být špatné. +more ** Mohou existovat lepší nelineární a/nebo zkreslené metody. ** Vlastnosti nejsou garantovány, pokud je použit nesprávný variogram. Nicméně obvykle je dosaženo ještě "dobré" interpolace. ** Nejlepší nemusí být nutně dobré: např. V případě neprostorové závislosti je krigovací interpolace jen tak dobrá jako aritmetický průměr. * Kriging používá \sigma_k^2 jako měřítko přesnosti. Nicméně toto opatření se opírá o správnost variogramu.
Aplikace
Ačkoli kriging byl původně vyvinut pro použití v geostatistice, je to obecná metoda statistické interpolace, která může být použita v jakékoliv disciplíně na vzorku dat z náhodných polí, které splňují příslušné matematické předpoklady.
K dnešnímu dni se kriging používá v celé řadě oborů, včetně následujících: * vědy o životním prostředí * hydrogeologie * těžba * přírodní zdroje * dálkový průzkum * hodnocení realit
a mnoho dalších.
Softwary, které využívají kriging
R packages # BACCO - Bayesian analysis of computer code software # tgp - Treed Gaussian processes # DiceDesign, DiceEval, DiceKriging, DiceOptim - metamodeling packages of the Dice Consortium.
* Matlab/GNU Octave # [url=http://mgstat. sourceforge. +morenet/]mGstat[/url] - Geostistics toolbox for Matlab. # [url=https://web. archive. org/web/20050604080848/http://www2. imm. dtu. dk/~hbn/dace/]DACE[/url] - Design and Analysis of Computer Experiments. A matlab kriging toolbox. # GPML - Gaussian Processes for Machine Learning. # [url=http://sourceforge. net/projects/kriging]STK[/url] - Small (Matlab/GNU Octave) Toolbox for Kriging for design and analysis of computer experiments. # scalaGAUSS - Matlab kriging toolbox with a focus on large datasets.
* Scilab # DACE-Scilab - Scilab port of the DACE kriging matlab toolbox # krigeage - Kriging toolbox for Scilab # KRISP - Kriging based regression and optimization package for Scilab
* Python # scikit-learn - machine learning in Python