Masarykova univerzita Přírodovědecká fakulta Geografický ústav Geostatistika doc. RNDr. Petr Dobrovolný, CSc. (pracovní texty) Brno 2010 i 1. Geostatistika - vymezení pojmu Geostatistika v užším slova smyslu - skupina interpolačních algoritmů založených na metodě krigingu. V širším slova smyslu - statistická analýza prostorově lokalizovaných dat. Pomocí „klasických" statistických metod lze vhodně analyzovat především atributová data - jejich kvantitativní či kvalitativní vlastnosti. Velmi omezeně však jimi lze charakterizovat prostorové vlastnosti objektů a jevů. Tyto prostorové vlastnosti jako např. spojitost jevů, prostorovou autokorelaci, prostorové uspořádání (strukturu) lze charakterizovat právě pomocí geostatistických metod. A B Cmínt 15251 15251 100 00 100 00 Standard Dcviation 20.00 20.00 \h 'Iííui 100.35 100 92 LO Pfín-eutile 7!!.mi 73.95 ÍW PwCTsntilR 125.61 124.72 Bíanple Dna Set A Direitcn: D.D Tolerance: &3.D zKinple Dala 3 CinaclÍDn: D.Q Tcleransa: 90.0 Obr. 1.1. Prezentace prostorového rozšíření spojitého jevu metodami popisné statistiky a pomocí tzv. semivariogramu. Na obrázku jsou znázorněny dva příklady zcela rozdílného prostorového rozšíření jistého spojitého jevu - např. koncentrace znečištění území jistou látkou. Z níže uvedené tabulky základních popisných 2 charakteristik i histogramů nelze zjistit žádný podstatný rozdíl v prostorovém uspořádání studovaného jevu v obou porovnávaných mapách. Ten je však patrný pokud prostorové rozšíření charakterizujeme pomocí tzv. semivariogramu, který patří k základním nástrojům strukturní analýzy a geostatistických metod. Geostatistika v širším slova smyslu představuje především: • Konstrukce spojitých polí tzv. deterministickými metodami • Koncept prostorové autokorelace • Strukturní analýzu a popis prostorové autokorelace strukturními funkcemi • Konstrukci spojitých polí metodami krigingu • Statistický popis prostorově lokalizovaných dat (geografických objektů) • Statistický popis prostorového uspořádání objektů (bodů, linií, ploch) • Objektivní metody klasifikace 2. Metody prostorové interpolace 2.1 Základní pojmy Prostorová interpolace slouží k odhadu hodnot určitého jevu či jeho intenzity v libovolném místě studované plochy, pro niž existují známé hodnoty tohoto jevu pouze v určitých lokalitách (meteorologické stanice, výškově zaměřené body apod.) Metod tedy lze využít ke konstrukci spojitých polí, k následné analýze prostorových dat - morfbmetrické a hydrologické modelování, optimální lokalizace apod.) Interpolace - skupina metod, které slouží k odhadu neznámých hodnot proměnné v jistých bodech (neměřených) na základě hodnot proměnné v bodech měřených. Prostorová interpolace - skupina metod, které slouží k vytváření spojitých povrchů (polí) z bodových měření. Body mohou být lokalizovány v 1, 2 i 3 rozměrném prostoru. Interpolace se může týkat nejenom bodů, ale i linií a ploch. V rámci interpolace je často řešen také problém extrapolace - tedy odhad hodnot proměnné vně oblasti definované krajními body měření. Naprostá většina interpolační ch postupů je založena na principu prostorové autokorelace - tedy na předpokladu, že hodnoty odhadované veličiny v lokalitách blízkých si boudou více podobné než hodnoty v lokalitách vzdálených. 2.1.1 Výběr reprezentativních vzorků Lokalizace měřených (odměrných) bodů v zájmovém území. Rozmístění (tzv. sampling) je důležité pro výběr interpolačního algoritmu a úspěšnost vlastní interpolace. Rozmístění • Pravidelné - může být zavádějící v případě zcela rovnoměrně rozmístěného jevu, který je studován (stromy, ...) • Náhodné - ze statistického hlediska je korektnější. Má ale i zápory - problematická lokalizace a zaměření jednotlivých míst než u pravidelně rozmístěných uzlů mřížky. Náhodné a nerovnoměrné rozmístění nemusí vystihovat základní rysy v rozložení měřené charakteristiky a může být i nákladnější. 3 a é t) rtgvttf t$ff!ptl!\g b) random sampling tí ŕrinsect ítmpling • • 1 • « •4 ■ • t • • \ • • • • • • • • • • \ • • • • c) strttifled rtmSem stmpUng cO cluster a*mptlag I) contour aimptlng Obr. 2.1 Možné způsoby rozmístění reprezentativních vzorků Jistým kompromisem mezi pravidelným a náhodným rozmístěním může být rozmístění stratifikované náhodné (stratified random). Shlukové uspořádání umožňuje studovat jev na několika měřítkových úrovních. V řadě případů je z různých důvodů (např. ekonomických, dostupnost, ...) prováděno měření pouze v omezené míře (profily - transepty). Ve velké části interpolačních úloh je rozmístění měřených bodů předem dáno, bez možnosti ho výrazně ovlivnit vhodnou lokalizací a výběrem (např. síť meteorologických stanic apod.) Prezentace spojitých polí - grid, TIN, izočáry, areály Možné datové zdroje pro interpolaci • bodová měření v terénu • digitalizované izolinie či polygony • stereopár leteckých fotografií či družicových obrazových záznamů Předpoklady úspěšné prostorové interpolace • existence dostatečně reprezentativního vzorku měřených dat • vhodné vlastnosti měřené veličiny a typ dat (ordinální, intervalová, poměrová) • teoretické i empirické znalosti o povaze prostorové diferenciace studovaného jevu • znalost podstaty použitelných interpolačních metod • znalost způsobu výběru nej vhodnej ší metody Běžné problémy interpolace: • vymezení studované plochy - přirozené a administrativní hranice • dostupnost bodů měření vně studované plochy 2.1.2 Průzkumová analýza dat (EDA - Exploratory Data Analysis) • ESDA - Exploratory Spatial Data Analysis • ESTDA - Exploratory Spatio - Temporal Data Analysis Množina statistických metod a speciálních nástrojů, zvláště grafických metod, používaných k lepšímu porozumění datům, k odhalení jejich důležitých vlastností. Jejím cílem je zjistit 4 základní informace o charakteru vstupních dat v tomto případě za účelem následné interpolace. Postupy a nástroje ESDA jsou využívány i v obecné prostorové analýze dat (studium prostorové autokorelace, pattern detectors). EDA slouží k průzkumu, deskripci, vizualizaci, zvýrazňování základních rysů dat, jejich distribuce (nejen ve smyslu prostorovém). Postupy EDA slouží k prověření požadavků normality, stacionarity vstupních dat. K těmto účelům používá specifických nástrojů (histogram, box plot, scatter plot, Q-Q graf). Deskriptívni metody používají jako měr úrovně ne „průměry" ale mediánu, počítají momenty vyššího řádu (asymetrie a špičatosti). Postupy EDA mohou vést k nutnosti úpravy či transformace původních dat. Úprava může spočívat v odstranění trendu či odlehlých hodnot, transformace potom např. například v tzv. log-transformaci. ESDA je nezbytným předstupněm úspěšné aplikace řady interpolační ch postupů (např. metod krigingu). Exploratory Spatial Data Analysis Selection of Data Points ._I_, Select by location Histogram tool Selection tool- ArcMap data view Voronoi mapping tool Select using ESDA tool ArcMap data view Histogram tool Ě-- i3Ja£ľ W Voronoi mapping too Obr. 2.2 Nástroje EDA jsou často propojeny s vlastní mapou (ESRI, Using ArcGIS Geostatistical analyst). Základní postupy průzkumové analýzy prostorových dat • výpočet základní popisné statistiky včetně momentů vyššího řádu (asymetrie a špičatosti) • prověření požadavků normality a stacionarity • analýza rozdělení hodnot - analýza histogramu • analýza kvantilového grafu (Q-Q grafu) • zkoumání odlehlých hodnot a jejich případné odstranění • analýza trendu a jeho případné odstranění • případná transformace vstupních dat (log) Základní nástroje ESDA Popisná statistika a „mapped histogram" - propojení mapy a grafu 5 aa at i or If 10 3D 30 « !(■ Odl <* **g vr* b*i S«t |« 3 f* Court :K SkMMi :*.*SH3 Ita ID Kurt™ 3 7« Mm -110 1 .i g«wti, <|>3 • Prediktor má podobu váženého průměru, tedy: + w23 + ... Doposud nebylo využito hodnot v měřených bodech. Proto váhy wi, W2 W3 jsou nalezeny na základě podmínky, že pokud je odhadován bod v bodě měření, je interpol ován přesně. Tato podmínka umožňuje sestavit soustavu N rovnic o N neznámých, která má jednoznačné řešení. Všechny metody interpolace využívající RBF dávají velmi podobné výsledky. Metody mají možnost nastavit parametr, který ovlivňuje shlazení výsledného povrchu. U všech metod RBF platí, že čím větší hodnota vyhlazovacího parametru, tím více shlazený je povrch. Opačně je tomu pouze pro tzv. inverzní multiquadric RBF. Nejčastěji využívanou je multikvadriková metoda (multiquadric RBF), která vychází z řešení následující rovnice: BI(x,y) = Vdi(x,y)2+R2 kde Bi(x,y) - radiální funkce vzdálenosti di(x,y) di(x,y) - relativní vzdálenost měřeného bodu v místě xi, yi od místa odhadu x, y R - vyhlazovací parametr Pro funkce Bi(x,y) jsou během výpočtu v každém interpol ováném bodě stanovovány váhy řešením soustavy lineárních rovnic. 2.4.8 Kriging Je to lokální interpolátor, který optimalizuje výběr bodů okolí, ze kterých je odhadována nová hodnota. K této optimalizaci se provádí tzv. strukturní analýza založená na studiu tzv. semivariogramu a konstrukci teoretického modelu. Parametry tohoto modelu jsou použity ve vlastním krigování. Kriging je založen na odhadu závislosti průměrné změny v hodnotách studované veličiny a vzdálenosti měřených bodů. Strukturní analýze a metodě korigování je věnována zvláštní kapitola. 2.4.9 Metody prostorové interpolace ploch (area based) Mnoho jevů se vztahuje k plošným jednotkám spíše než k bodům (hustota obyvatelstva států, kvalita pitné vody...). Metody řeší způsob, jakým lze odhadnout hodnoty jistého jevu na základě hodnot jiného jevu vázaných na plošné jednotky. Přitom mohou nastat dvě situace: 1. plošné jednotky se shodují 2. zdroj ové j ednotky j sou podmnožinou (nested) j ednotek výstupních Metody lze dělit do dvou skupin: 27 1. metody zachovávající objem studovaného jevu (volume preserving) 2. metody nezachovávající objem studovaného jevu (non-volume preserving) Metody nezachovávající objem studovaného jevu (non-volume preserving) Příklad - mapa A vyjadřuje celkový počet obyvatel ve čtyřech administrativních jednotkách určitého území. Mapa B potom záplavovou oblast kolem vodního toku. Cílem je zjistit, jaká je hustota obyvatelstva uvnitř záplavové zóny. Map*. Administrative boundaries with Map H: Hood hazard map tor total Martian population counts it? canal A B 200 400 □ C 600 eoo canal flood zona Map C: Population Area oacti call im *2 "4 *6 'a Map D: Interpolated population density grid 2 3 3 4 5 S 4 5 7 6 • e 7 7 Map E; Einmal td Mil ptyiMrm for grid SO 75 100 100 75 100 125 125 100 125 175 150 160 150 175 j 175 Map F: Total Martian population counts lor Hood zone Obr. 2.28 Princip metody prostorové interpolace ploch nezachovávající objem studovaného jevu Postup: 1) výpočet hustoty obyvatelstva pro každou plochu 2) určení centroidu každé plochy 3) interpolace hustoty obyvatelstva výše popsanými metodami Metody zachovávající objem studovaného jevu (volume preserving) Provede se překrytí cílových zón (oblastí) přes oblasti zdrojové a určí se poměrná část cílové zóny, která spadá do zóny zdrojové. Celková hodnota atributu v cílové zóně je určena v závislosti na plošném zastoupení zón zdrojových. 28 2.4.10 Pycnophylatic method Jeden z hlavních problémů metody thiessenových polygonů je, že měřený prvek se považuje za homogenní v rámci jedné třídy, veškeré prostorové změny jsou vázány na hranice. V případě spojitých či relativně spojitých prvků jde o naprosto nevhodný způsob interpolace. Modifikací je metoda, která v rámci každé třídy zachovává sumu studovaného prvku, avšak dovoluje kontinuální změnu směrem k hranicím každé třídy. Metoda bere v úvahu hodnotu atributu sousedních tříd. Nepředpokládá se existence bariér a hodnoty sousedních tříd jsou shlazeny bez skokových změn v hodnotách daného atributu pomocí pravděpodobnostní funkce (density function). Metoda byla použita na demografická data. Jde o neexaktní interpolator. Minimální i maximální hodnoty vypočtené touto metodou j sou výrazně menší resp. větší než skutečně naměřené hodnoty. Map A: Administrativ* bound! liwwith total Martian population: counts A B 200 400 D C 60O BOO M*p B: Fl»d hazgnj m?p far the cans. canal flood zone Map C; Area ratio population ccunis 100 200 100 200 300 400 300 405 j D: Estimated total population counts fot flood zona 300 1000 flood -haZaM ZOflti TOO Obr. 2.29 Princip metody prostorové interpolace ploch zachovávající objem studovaného jevu Radiál basis interpolation - skupina metod podobných geostatistickým metodám krigování, nrlze ale využít metody modelování variogramu, nejsou vyžadovány vstupní podmínky jako u korigování. Splines SPLÍNE Cubic - polynom 3 řádu 29 Bicubic - 2 souřadnice ř(x,y)=ai +a2x+a3y+a4x2+a5ý'+a6xy+a7x2y+ae,xyq+a9/+a^]ř Remember that linear interpolation uses a linear function for each of intervals [xk,Xk+i]. Spline interpolation uses low-degree polynomials in each of the intervals, and chooses the polynomial pieces such that they fit smoothly together. The resulting function is called a spline. For instance, the natural cubic spline is piecewise cubic and twice continuously differentiable. Furthermore, its second derivative is zero at the end points. The natural cubic spline interpolating the points in the table above is given by 3. Geostatistické metody interpolace Při použití většiny dosud uvedených interpolační ch algoritmů nemáme a priori objektivní informaci o tom, zda způsoby definování okolí interpolovaného bodu, body použité k interpolaci a také jejich váhy jsou zvoleny vhodně. V případě kvalitních vstupních dat dává většina interpolační ch technik podobné výsledky (Borrough, McDonnell, 1998). Vstupních dat musí být především dostatečný počet a musí být také vhodně rozmístěna ve studovaném území. V opačném případě hraje velkou roli volba vhodného interpolačního algoritmu. Žádná z metod dosud neřešila následující problémy: • počet bodů nutných k výpočtu lokálního průměru • velikost orientaci a tvar okolí • zda neexistuje jiná cesta k definování vah než funkce vzdálenosti bodů • jaké jsou chyby a nejistoty spojené s interpolovanými hodnotami Odpovědi na tyto otázky poskytují geostatistické postupy založené na tzv. strukturní analýze. Její výsledky jsou poté mimo jiné využitelné v interpolačních postupech krigingu. Metodu krigingu uvedli G. Matheron a D.G. Krige a je založena na skutečnosti, že interpolovaný povrch lze lépe vystihnout stochastickou funkcí než shlazující matematickou funkcí. Předchozí metody interpolace lze označit jako deterministické. Metoda krigingu odhaduje interpolační váhy bodů tak, že optimalizuje interpolační funkci a lze ji například využít také k optimalizaci sítí apod. Kriging jako interpolační metoda byla vyvinuta v geologii a je vhodná pro interpolovaní proměnných, které se v prostoru mění s jistou kontinuitou, ale nelze je popsat jednoduchou shlazující funkcí některého z globálních interpolátorů. V případě krigingu je interpolovaný povrch tvořen ze tří složek (obr. 1). První složku představuje tzv. strukturální komponenta s konstantním průměrem či trendem - obecný trend (tzv. drift). Druhou složku povrchu představují kolísání (drobné sní ženiny či vyvýšeniny), jejichž podstatu lze vyjádřit určitou matematickou funkcí jako v případě trendu, ale která vyjadřují určitou prostorovou korelaci (tzv. regionalizovaná proměnná) Třetí složku potom představuje náhodný šum. 30 Obr. 3.1 Základní komponenty spojitého povrchu (i - trendová složka - drift; ii - tzv. regionalizovaná proměnná; iii - náhodná složka) V případě krigingu jsou všechny tři složky analyzovány separované. První složka je odhadována za pomoci obecné trendové funkce. Druhá složka - náhodná, ale prostorově korelovaná kolísání jsou analyzována metodou tzv. variogramu. Nechť x je poloha bodu v 1, 2 či 3 dimenzích, potom hodnota náhodné proměnné Z v bodě x bude z(x) = //(x) + s' (x) + s x - pozice v 1, 2 či 3 rozměrném prostoru Z - interpolovaná proměnná Z(x) - hodnota proměnné v bodě x ju(x) - deterministická složka (trend) s'(x) - stochastická složka (regionalizovaná proměnná) - lokálně proměnné, ale prostorově závislé reziduum od ju(x) s" - náhodná, prostorově nezávislá složka, gaussovský šum s nulovým průměrem a s rozptylem o . Velké písmeno Z značí, že se jedná o náhodnou funkci a ne o měřenou hodnotu proměnné z. Prvním krokem je zvolení vhodné funkce vystihující složku ju(x). V nejjednodušším případě, když se v hodnotách nenachází žádný trend (drift), potom se fi(x) rovná průměrné hodnotě v ploše a průměr či očekávaný rozdíl mezi dvěma místy x a x+h vzdálenými od sebe o hodnotu vektoru h, bude nula: E[Z(x)-Z(x+h)] = 0 kde Z(x) a Z(X+h) jsou hodnoty náhodné proměnné Z v poloze x, x+h. Dále také předpokládáme, že rozptyl rozdílů závisí pouze na vzdálenosti mezi místy, tedy: e[{Z(x) - Z(x + h)}2 \ = e[{s' (x) -s\x + h)f\= 2y(h) kde hodnota y(h) se označuje jako semivariance. 31 To značí, že jakmile jsme odhadli příspěvek strukturní proměnné ju(x), zbývající kolísání má konstantní rozptyl a diference mezi dvěma místy jsou pouze funkcí jejich vzdálenosti. Výše uvedený vztah lze přepsat: Z(x) = //(x) + r(h) + e" tedy mezi s '(x) a y(h) je ekvivalence 3.1 Strukturní analýza Geostatistická strukturní analýza (variografie) je procedura zahrnující výpočet strukturálních funkcí, výběr a konstrukci odpovídajících teoretických modelů a jejich aplikace, interpretaci průběhu strukturálních funkcí. Cílem je popsat takové vlastnosti jako jsou kontinuita, homogenita, stacionarita či anizotropie pole studovaných prostorových proměnných veličin. Tyto vlastnosti jsou popisovány prostřednictvím měr prostorové autokorelace a prostorové variability. Strukturální analýza je výchozím krokem geostatistického modelování Sama o sobě ale poskytuje řadu velmi důležitých informací o struktuře náhodného pole jako modelu konkrétního objektu v krajinné sféře. Obr. 3.2 Příklad výpočtu měr prostorové variability pro ID (řadu hodnot) Ke kvantifikaci prostorové autokorelace, která vyjadřuje skutečnost, že objekty blízké si jsou více podobné než objekty vzdálenější slouží strukturální funkce. Strukturální funkce měří sílu korelačního vztahu jako funkci vzdálenosti. Na obr. 2 je pro j ednoduchost odvozena strukturální funkce pro řadu pravidelně rozmístěných bodů na linii (tedy pro ID). Charakteristiky, které popisují prostorovou variabilitu lze odvodit z měr úrovně a variability následovně: průměr = (1 +3+6+5+3+1+2+3)/8=3,0 rozptyl= [(1 -3)2+ (3-3)2+ (6-3)2+ (5-3)2+ (3-3)2+ (1 -3)2+ (2-3)2+ (2-3)2] /8= 2,75 kovariance(l)=[(l-3)*(3-3)+(3-3)*(6-3)+(6-3)*(5-3)+(5-3)*(3-3)+(3-3)*(l-3)+(l-3)*(2-3)+(2-3)*(3- 3)]/7=l,14 semivariance(l )= [ (1 -3)2+ (3-6)2+ (6-5)2+ (5-3)2+ (3-1)2+ (1 -2)2+ (2-3)2] 11= 3,43 semivariance(2)= [(l-6)2+ (3-5)2+ (6-3)2+ (5-l)2+ (3-2)2+ (l-3)2]/6=9,83 semivariance(3)=[(l-5)2+(3-3)2+(6-l)2+(5-2)2+(3-3)2]/5= 12,50 Semivariance může být definována jako. 32 /(xixj) = l/2var(Z(xi)-Z(xj)) kde var značí rozptyl. Jestliže jsou dva body xi a Xj blízko sebe, bude rozdíl hodnot studované veličiny Z(xO a Z(xj) těchto bodech malý. S růstem vzdálenosti si budou hodnoty méně podobné. Grafickým vyjádřením závislosti semivariance na vzdálenosti je strukturální funkce nazývaná semivariogram, jejíž typický průběh je na obr. 3. Obr. 3.3 Vztah mezi semivariogramem (1) a kovarianční funkcí (2). Vysvětlení jednotlivých termínů viz. dále Semivariogram je pouze jednou z mnoha strukturálních funkcí, i když nej používanější. Semivariogram je mírou nepodobnosti. K dalším takovýmto funkcím patří kovarianční funkce: C(xixJ) = cov(Z(xi),Z(xJ)) kde cov značí kovarianci. Kovariance je stejně jako korelace mírou podobnosti. Budou-li dva body blízko sebe, budou hodnoty v těchto bodech podobné a jejich kovariance bude velká. S růstem vzdálenosti bude klesat podobnost bodů a budou klesat i hodnoty kovariance. Hodnota kovariance je při nulové vzdálenosti rovna rozptylu množiny zpracovávaných dat. Vztah mezi semivariogramem a kovarianční funkcí lze vyjádřit následovně: /(xixj) = siU-C(x;,xj) Veličina označená sill značí tzv. práh a je to hodnota, na níž má semivariogram vodorovný průběh (viz. dále). 3.2 Strukturní analýza v2D Výše uvedený příklad určení semivariance platí pro řadu pravidelně rozmístěných bodů. Strukturní analýzy a následného krigování se však ve většině případů používá pro charakterizování vztahů prostorové autokorelace a pro následné odhady a interpolace v rovině, ve kterém máme množinu pravidelně, častěji však nepravidelně rozmístěných bodů měření. Abychom popsali závislost hodnot studovaného jevu v prostoru, vyneseme hodnoty semivariancí pro všechny dvojice bodů do semivariogramu obdobně jako ve výše uvedeném případě. Semivariogram je strukturální funkce, která popisuje závislost průměrné kvadratické diference hodnot prostorové proměnné veličiny Z na vzdálenosti h. Semivarianci lze odhadnout z naměřených dat podle následujícího vztahu: 33 rth^fžízoo-zfc+h))2 kde: n - počet dvojic bodů pozorování proměnné s atributem z vzdálených o hodnotu h h - tzv. lag - vzdálenost dané dvojice bodů. h Obr. 3.4 Experimentální semivariogram (+) s charakteristickými hodnotami pro vzdálenosti h (•) a proložený teoretický model semivariogramu (plná čára) Graf hodnot /(h) a h se označuje jako experimentální semivariogram a je prvním krokem ke kvantitativní deskripci regionalizované proměnné. I když je semivariogram vektorová funkce, sestavuje se často jako izotropní (tj. bez ohledu na směr) pro celkové charakterizování hodnoceného náhodného pole nebo tehdy, je-li k dispozici omezený počet pozorování. Experimentálním semivariogramem se v následujícím kroku prokládá teoretický model. 3.2.1 Prvky semivariogramu Na Obr. 3.5 je uveden často používaný tzv. sférický model s vysvětlením používané terminologie. Na horizontální oseje vynášena vzdálenost h mezi jednotlivými vstupními body interpolovaného povrchu (tzv. lag), na vertikální ose potom rozptyl zkoumané proměnné jako funkce vzájemných vzdáleností jednotlivých měřených bodů. Je mírou, která vyjadřuje, jak velké je okolí, daného bodu, ve kterém se nacházejí body sousední, jejichž hodnota interpolovaného atributu závisí (koreluje) s hodnotou v tomto bodě. Takto vynesenými body je proložena křivka mající charakteristický tvar. Je-li vzdálenost mezi dvěma body malá, jejich hodnoty jsou podobné a hodnota semivariance je také malá. Se zvětšující se vzdáleností hodnota semivariance roste. Při určité vzdálenosti dvou boduje možno říci, že jejich hodnoty (např. výšky) spolu již nekorelují a hodnota semivariance i se zvyšováním vzdálenosti již neroste, ale zůstává konstantní. Plochá část semivariogramu určuje tzv. práh (sill) a je rovna rozptylu zpracovávaných dat. Existence prahu značí, že se zvětšující se vzdáleností se hodnoty semivariance nemění. Kritická hodnota vzdálenosti, na níž se křivka semivariogramu stává rovnoběžnou s vodorovnou osou se označuje jako dosah (range). Dosah definuje pro daný bod velikost okolí, které je nutné uvažovat při interpolaci hodnoty v daném bodě. 34 d i ✓ C1 t/ s2 C:, Obr. 3.5 Příklad teoretického semivariogramu - sférický model. Parametry semivariogramu: a - dosah (range), d - rozpětí, c0 - zbytkový rozptyl (nugget), c=c0 + ct - práh (sill), h - lag (krok vzdálenosti) Z výše uvedeného by také mělo platit, že proložená křivka semivariogramu by měla procházet počátkem souřadné soustavy (nulová vzdálenost mezi body znamená zákonitě také nulovou hodnotu rozptylu). Velmi často nenabývají experimentální semivariogramy v počátku nulové hodnoty; protínají osu y v nenulové hodnotě, která je nazývána zbytkový rozptyl (nugget variance). To může ukazovat na rozptyl menší než je "vzorkovací" vzdálenost, nebo na malou přesnost měření, kdy např. jsou v datech obsaženy dva vzorky ze stejného místa, pokaždé s jinou hodnotou. Zbytkový rozptyl je tak vyjádřením náhodného šumu s" a sestává se jednak z chyb měření (dvě různé hodnoty projeden bod) ajednak z tzv. microscale variation - ty jsou vyjádřením rozptylu hodnot složky s'. V případě nulového zbytkového rozptylu a tedy v případě nulové chyby měření je krigování exaktním interpolátorem. Na tvar semivariogramu má značný vliv velikost hodnoty lag (vzdálenosti h). Velikost h se volí např. jako průměrná minimální vzdálenost mezi sousedními body. Velká hodnota h dává hladší průběh semivariogramu, může však zamaskovat efekt autokorelace studovaných dat na menší vzdálenosti. Naopak malá hodnota h má vliv na malý a tedy často nereprezentativní počet bodů v rámci každé hodnoty násobku h. Nevhodná délka kroku se může dále projevit v oscilaci hodnot semivariogramu. Výpočet semivariancí se většinou provádí do vzdálenosti rovné polovině maximální vzdálenosti bodů v prostoru. Tedy násobíme-li velikost kroku počtem kroků, měli bychom dostat hodnotu rovnou zhruba polovině maximální vzdálenosti mezi interpolovanými body. 3.2.2 Efekt anizotropie V případě velkého množství nepravidelně rozmístěných boduje vhodné hodnotu semivariance vyjádřit pro skupiny bodů přibližně stejně vzdálených. V tomto případě se do výpočtu hodnoty semivariance pro dané h berou všechny body padnoucí do mezikruží určeného tolerancí délky kroku. Toleranci délky krokuje nutné volit v případě nerovnoměrného rozmístění měřených bodů v rovině. Tato hodnota se volí v mezích od 10 -50% z délky kroku h. Grupovaní hodnot semivariancí na základě podobné vzdálenosti (tzv. binning) dovoluj e konstruovat druhý typ grafu, často využívaného pro studium prostorové autokorelace a pro snadnější interpretaci hodnot semivariogramu - tzv. plošný graf semivariance. Ten navíc umožňuje posoudit eventuelní rozdíly v hodnotách semivariance 35 v závislosti na směni - tedy definovat efekt anizotropie. Proto se tento typ grafu také označuje jako povrch anizotropie. Plošný graf semivariance představuje grid s buňkami o velikosti strany rovné vzdálenosti h (lag). V grafu určujeme vzdálenost směrem od středu. Jednotlivé buňky grafu nesou hodnotu semivariance vypočtenou ze všech dvojic bodů, které jsou od sebe vzdáleny o právě o vzdálenost od středu grafu a které se navíc nacházejí v určitém směru - viz. obr 6. 2) a 3). Hodnoty semivariancí jsou potom vyjádřeny např. barvou. 2 Obr. 3.6 Princip grupovaní hodnot semivariancí na základě podobné vzdálenosti a plošný graf semivariance. Hodnoty semivariancí obecně rostou směrem od středu grafu, protože podobnost hodnot studované veličiny s růstem vzdálenosti obecně klesá - tedy roste jejich nepodobnost vyjádřená semivariancí. V případě, že se hodnota semivariance mění s rostoucí vzdáleností (směrem od středu grafu) stejně ve všech směrech, potom hovoříme izotropii studovaného pole. V opačném případě tvoří hodnoty semivariance na plošném grafu tvar elipsy. Výsledkem je, že i tvar semivariogramu bude jiný ve směru hlavní a vedlejší poloosy. Semivariogram sestavený z bodů ve směru kratší poloosy elipsy anizotropie se bude vyznačovat strmějším průběhem. Tento směr odpovídá směru maximální variabilizy manimálníhu dosahu (viz. dále). Směr hlavní osy elipsy je naopak směrem minimální variability. 36 Obr. 3.7 Povrch vykazující efekt anizotropie a odpovídající empirické semivariogramy Tzv. izotropní semivariogram tedy neuvažuje odchylky v závislosti na směru, naproti tomu anizotropní semivariogram se liší především odlišnou hodnotou dosahu pro specifické směry, další charakteristiky semivariogramu (typ, práh, zbytkový rozptyl) se většinou nemění. Takovouto anizotropii označujeme jako geometrickou. V případě, že nelze použít stejný model semivariogramu resp. stejné hodnoty prahu a zbytkového rozptylu hovoříme o tzv. zonální anizotropii. K modelování zonální anizotropie lze využít konstrukce tzv. složených modelů semivariogramu. Ke konstrukci směrových semivariogramu je nutné řešit otázku vhodného výběru bodů. Ve většině případů se volí 4 až 8 směrů a k jejich vymezení je nutné stanovit následující parametry (jejich význam je zřejmý z následujícího obrázku): • úhlovou toleranci • šířku pásma • délkovou toleranci (lag) šířka pásma Obr. 3.8 Parametry tzv. směrových semivariogramu 37 Efekt anizotropie je vyjádřením náhodného procesu chovaní studované veličiny. Nelze ho zaměňovat s trendovou složkou. Ta by měla být ve zpracovávaných datech předem definována a před vlastní strukturní analýzou odstraněna (viz. ESDA). Pokud k tomu není pádný důvod, daný fyzikální podstatou zpracovávaných dat, není vhodné používat anizotropního modelu s poměrem os elipsy anizotropie větším než 3 ku 1. Pokud experimentální semivariogram ukazuje na takovou anizotropii, je to způsobeno zřejmě trendem obsaženým v datech. Potom je vhodné nejprve z dat trend odstranit. Teoretický semivariogram Jedná se o model, který nejlépe aproximuje průběh experimentálního semivariogramu v okolí počátku a prahu (viz. dále). Právě proces hledání teoretického semivariagramu se někdy označuje jako strukturní analýza. Modely semivariogramu se dělí podle chování v okolí počátku a v „nekonečnu" do několika skupin: • modely přechodového typu - tj. s prahem (sférický, kvadratický, gaussovský, exponenciální), • modely bez přechodu (lineární, logaritmický), • modely s oscilujícím prahem (sinový, cosinový), • čistě náhodný model. U prvních tří skupin se může objevit tzv. efekt zbytkového rozptylu (nugetový efekt), který se odráží v posunu grafu semivariogramu o hodnotu c0 ve směru osy /(h). Model nalezený pro danou množinu dat závisí jak na experimentálních, tak teoretických předpokladech. Vlastnosti, které prakticky vedou k určení konkrétního teoretického modelu, jsou: • přítomnost nebo absence "ploché části" semivariogramu - tzv. prahu; v rovnicích semivariogramu je dán konstantou C • vzdálenost, ve které semivariance dosáhne prahové hodnoty - dosah (range); v rovnicích semivariogramu je dán konstantou a • chováním v počátku (tj. semivariance mezi velmi blízkými body) • dosah je mírou korelace uvnitř množiny dat; "dlouhý" dosah indikuje vysokou korelaci, "krátký" dosah korelaci nízkou. • hodnota prahu je rovna celkovému rozptylu. 3.3 Přehled teoretických semivariogram 3.3.1 Modely přechodového typu Modely přechodového typu (transitivní) - prostorová autokorelace kolísá s hodnotu h. U těchto klasických modeluje vyjádřena skutečnost, že při malých vzdálenostech je shoda mezi zjištěnými hodnotami vysoká (a tedy variabilita nízká), s rostoucí vzdáleností se „neshoda" zvyšuje až do určité vzdálenosti (dosah), kde se úroveň neshody stabilizuje kolem hodnoty 38 statistického rozptylu. Za touto vzdáleností se již neuplatňuje prostorová vazba mezi zkoumanými místy a variabilita je plně určována statistickým rozptylem. Sférický model - zbytkový rozptyl je důležitý, ale malý. Je jasně vyjádřen dosah (range) a prahová hodnota (sill). Je typický pro pole, ve kterém dominuje jeden zdroj variability. d s2 Obr. 3.9 Sférický model semivariogramu >3~ r(h) = c0 + ct * 3h_ 1 h 2 a 2I a ........... pro h < a y(h) = c0 + Cj ...........pro h > a Kvadratický model d 1 i s2 yk(») Co a = 2d Obr. 3.10 Kvadratický model semivariogramu y(h) = c0 + c, * "2i-a '0 Cl pro h < a Exponenciální model - dobře vyjádřené hodnoty zbytkového rozptylu a prahu, ale pouze postupná aproximace k hodnotě dosahu (range) 39 y*(h) d [- S2 y y*(« Co a = 3d Obr. 3.11 Exponenciální model semivariogramu y(h) = c0 + q * [l-exp(-h/d)], kde a = 3d Gaussův model - hladký povrch, hodnota zbytkového rozptyluje velmi malá ve srovnání s regionalizovanou proměnnou. Model má inflexní bod. Je typický plynulými změnami hodnot. Používá se často např. při modelování výškových dat. Je používán u dobře prozkoumaných polí. Často se však vyznačuje nestabilitou. y*(h) d ¥ Co a h Obr. 3.12Gaussův model semivariogramu y(h) = c0 + c * [l - exp(-h2 ld2)\, kde a = dyfŠ Lineární model s prahem - jednoduchý a poměrně často využívaný zvláště programy provádějícími interpolaci pomocí krigování na základě automaticky vypočítaného a vyhodnoceného semivariogramu. Při provádění strukturální analýzy se využívá raději jiných přechodových modelů. y*(h) d Co a h Obr. 3.13 Lineární model semivariogramus prahem y(h) = c0 + bh ...........pro h < a 40 y (h) =c0 + c1 ...........pro h > a 3.3.2 Modely bez přechodu Modely bez přechodu (netransitivní) - nemají prahovou hodnotu (sill) v rámci studované plochy a lze je popsat např. lineárním modelem. Výskyt těchto modelu si lze zjednodušeně představit jako určitý extrémní případ klasického přechodového modelu. Představme si, že bychom u něho prováděli výpočet semivariogramu jen do vzdálenosti nepřesahující rozpětí d. Pak bychom při vynesení bodů nenašli žádnou oblast stabilizace křivky semivariogramu a daný případ bychom interpretovali jako model bez přechodu. Lineární model: h, h2 Obr. 3.14 Lineární model semivariogramu y(h) = c0 + bh, kde b je směrnice přímky Logaritmický semivariogram ln(h,) ln(h,) Obr. 3.15 Logaritmický model semivariogramu y(h) = c0 +bln(h), kde b je směrnice přímky 3.3.3 Oscilační modely Oscilační modely - oscilační (tj. nehomogenní) charakter má zkoumané pole nejčastěji v důsledku pravidelného střídání pásů s vyššími a nižšími hodnotami. Průměrná šířka pásů se dá odhadnout podle rozměru poloviny periody vlny. U těchto modelů se často projevuje nestabilita. Nepoužívají se pro odvození parametrů potřebných pro krigování (upřednostňují se robustní, jednoduché přechodové modely). 41 Jev, kdy hodnoty semivariagramu v jisté vzdálenosti delší než dosah začnou opět klesat či vykazují více lokálních minim ukazuje na periodická kolísání v hodnotách atributu a označuje se jako hole effect. Sinový model semivariogramu d 1 —\ s3 Co h Obr. 3.16 Sinový model semivariogramu sin(gh) r(h) = c0 + c, 1 gh kde g = nl co Hodnota sin se udává v radiánech. Dochází k postupnému tlumení hodnot oscilací. Hodnota co udává průměrný rozměr bohatších a chudších úseků. Cosinový model semivariogramu - Nedochází k postupnému tlumení hodnot oscilací. rm Obr. 3.17 Cosinový model semivariogramu y(h) = c0 + Cj[l - cos(gh)] kde g=7ľlco Čistě náhodný model semivariogramu r(h) = c( Semivariogram nemá žádnou úvodní rostoucí větev, hodnoty často pouze kolísají kolem prahu. K této situaci dochází, když je studované pole příliš variabilní vzhledem ke zvolenému kroku vzorkování (zjišťování hodnot). 42 3.4 Další druhy semivariogramů 3.4.1 Složené modely (komplexní semivariogram) 0,18 n 0.16 - g 0.12 - .§ 0.1 - .i 0.08- i 0.06 - 0 0 2 4 ô 10 Lag (m) Obr. 3.18 Složený model semivariogramů rT(h) = ^(h) + r2(h) + r3(h) +... Prostorová kolísání v závislosti na odlišných typech povrchů (cover classes) - svoji vlastní strukturu prostorového uspořádání a autokorelace hodnot proměnné mohou mít rozdílné kategorie landuse, druhy půd, atd. V tomto případě mohu být modely sestavené pro jednotlivé třídy vhodnější než model globální. Je zde však často problém dostatku dat. Indikátorové semivariogramy se konstruují a využívají při strukturální analýze nominálních (kvalitativních) dat (barva, druh horniny). Primární data se transformují do hodnot 1 a 0 podle splnění indikační podmínky - např. zdaje hornina pískovcem. Často slouží jako vstup pro tzv. indikátorové krigování (viz. dále). Soft semivariogramy se využívají při v případě nedostatku primárních dat, kdy je možné na základě provedené simulace doplnit další data a usnadnit provedení strukturální analýzy. Interpretace a verifikace je však dosti nesnadná a vyžaduje větší zkušenosti. Soft semivariogramy se často používají při provádění soft krigingu (viz. dále). 3.5 Analýza a interpretace strukturálních funkcí Pro každý model existují vlastní pravidla interpretace. Konstrukci semivarioagramu a odvození teoretického modelu by měla vždy předcházet důkladná analýza vstupních dat založená na metodách popisné statistiky (ESDA - explorační analýza prostorových dat, viz. Pro korektní odhady vhodného teoretického modeluje důležitý počet bodů uvažovaných pro vyjádření hodnot semivariance pro daný lag (h). Proto se často hodnoty teoretického modelu odhadují za pomoci vážené metody nejmenších čtverců, kdy jako váhy se berou počty párů na dané vzdálenosti h. Značný podíl šumu ve variogramu může být dále způsoben malým rozsahem vzorku použitého k výpočtu /(h) dále) 43 K dosažení stabilních hodnot se doporučuje 20 - 30, v některých případech však až až 50-100 hodnot. Je-li jejich počet nízký, stoupá chyba odhadu. Hladší průběh semivariogramu lze docílit zvětšením velikosti vyhledávacího okna (větším h). O velikosti okna vypovídá hodnota dosahu (range). Je-li odhadnutý dosah z variogramu příliš malý a všechny body jsou dále jak dosah, potom nejlepším odhadem je použití celkového průměru. Vzdálenosti mohu být modifikovány efektem anizotropie - potom je nutné měnit tvar okolí. Anizotoropie však může být výsledkem i nedostatečného počtu vzorků. Výpočet experimentálních semivariogramu se doporučuje provádět do vzdálenosti h < L/2, kde L je maximální vzdálenost míst pozorování v poli. Vždy je vhodné upřednostňovat jednodušší teoretický model semivariogramu, který dobře vystihuje hlavní rysy experimentálních hodnot, před modelem složitějším. V případě výpočtu experimentálního semivariogramu z nepravidelně sítě pozorování je nutno počítat s vyšší „rozkolísanosti" stanovených bodů kolem teoretického modelu. Úroveň prahu se obvykle doporučuje volit podle hodnoty statistického rozptylu. Je-li hodnota dosahu použitého teoretického semivariagramu malá vzhledem k hodnotám empirickým je možné zmenšit hodnotu kroku h a naopak Při prokládání tečny počátkem experimentálního semivariogramu pro určení rozpětí musíme respektovat skutečnost, že funkce semivariogramu je vždy kladná. Hodnota rozpětí je důležitá pro aplikaci oscilačních semivarigramů. Při interpretaci zbytkového rozptylu musíme uvážit i možný vliv chyb měření (technických chyb) výchozích pozorování. Výběr vhodného teoretického modelu musí vycházet z cíle analýzy. Je-li cílem odhalení strukturálních úrovní a podrobný popis všech charakteristik studovaného pole, pak je nutno podrobně analyzovat chování v celém reálném průběhu experimentálního semivariogramu. Jestliže je interpretace prováděna pro účely návazných krigovacích výpočtů, je účelné zvolit pokud možno jednoduchý a robustní model, vystihující chování a okolí počátku až do úrovně prahu. Při interpretaci j e důležité vycházet z dobré znalosti objektu v krajinné sféře a z využití všech informací o jeho parametrech. Při analýze anizotropie je podle zkušenosti dobré volit pro všechny směrové semivariogramy - samozřejmě pokud je to možné - stejný teoretický model. Proto je výhodné vyjít z izotropního semivariogramu pole. V případě anizotropního pole se zpravidla snažíme využít předpokladu geometrické anizotropie, kterou lze snadno eliminovat transformací souřadného systému. 44 Obecně je účelné postupovat tak, že v počáteční fázi aplikace geostatistických metod na přírodní objekt se provede podrobná interpretace strukturálních funkcí a v následných fázích se podle získaných zkušeností použije zjednodušený základní model. Analýza semivariogramu je podstatným krokem k určení optimálních vah pro interpolaci. Jestliže ve semivariogamu dominuje náhodná složka (s"), potom data obsahují takový šum, že interpolace nemá smysl. Jako nejlepší odhad z(x) je vhodné použít průměrnou hodnotu. Charakteristiky pole popsané strukturní analýzou: Kontinuita - je vyjádřena hodnotou dosahu semivariogramu. Pole s větší kontinuitou se vyznačuje vyšší prostorovou autokorelací. Nehomogenita - projevuje se tzv. oscilací hodnoty prahu. Délka poloviny periody odpovídá průměrnému rozměru elementů nehomogenity. Nehomogenity na dané úrovni pozorování nepostižitelné se projeví jako zbytkový rozptyl. Nestacionarita - projevuje se zpravidla parabolickým nárůstem křivky semivariogramu. Prokazatelná je případech, kdy dochází k parabolickému růstu křivky až za hodnotou dosahu, tedy na stabilizované části křivky. Nestacionarita pole dokládá změnu průměrné hodnoty proměnné v poli. Ze vzdálenosti, kde se začne deformace křivky semivariogramu projevovat, lze určit vzdálenost, do které jsou změny průměrné hodnoty v poli zanedbatelné. Anizotropie - lze ji popsat pomocí modelů jednotlivých směrových semivariogramu (tj. semivariogramu vypočtených na různých směrech v poli). Projevuje se změnami parametrů (dosahu, prahu, zbytkového rozptylu), jednak v rozdílech typů směrových semivariogramu. Jak bylo uvedeno výše rozlišujeme geometrickou a zonální anizotropii (viz. obr). Obr. 3.19 Rozdíl mezi geometrickou (A) a zonální (B) anizotropii semivariogramu 4. KRIGING - geostatistické metody interpolace Krigování je základní geostatistickou metodou určování lokálního odhadu. Metoda se často označuje akronymem BLUE (Best Linear Unbiased Estimator - tedy nejlepší lineární nezkreslený odhad). Toto označení má vystihnout výchozí podmínky krigování: • odhadovaná hodnota j e vypočtena j ako lineární kombinace vstupních hodnot: 45 n ž(xo)=Z!/íi-z(xi) i=l kde pro váhy platí í=i • nezkreslený (nestranný) odhad značí, že průměrná chyba tohoto odhadu je rovna nule • je minimalizován rozptyl odhadu ^(žj -Zj)2 = min. Pokud prostorově závislá náhodná kolísání nejsou překryta nekorelovaným šumem, potom může být semivariogram využit k určení vah X{ potřebných pro interpolaci. Procedura je podobná jako v případě metody vážených klouzavých průměrů s tím rozdílem, že právě váhy jsou odhadnuty geostatistickými metodami. Váhy jsou zvoleny tak, aby odhad ž(x0) byl nestranný a odhad rozptylu je tzv. Lagrangeův multiplikátor, který zajišťuje požadavek minimalizace odchylek a zároveň podmínku, že suma vah je rovna jedné. Uvedená metoda se označuje jako základní (ordinary) kriging a je možněji použít pro interpolaci v pravidelné mřížce hodnot, ke konstrukci map (např. izolinií). 46 PŘÍKLAD: Výpočet neznámé hodnoty v bodě metodou základního krigingu Na základě změřených hodnot veličiny Z v pěti bodech (i=l,..., 5) (viz. obrázek) máme za úkol odhadnout hodnotu Z bodě (i=0) o souřadnicích (x=5, y=5) metodou krigingu. Obr. 4.1 Vstupní data pro lokální odhad metodou základního krigování (podtržená čísla značí hodnotu atributu v bodě) Na základě předem provedené strukturní analýzy použijeme sférický semivariogram. y (h) = c0 + Cj * 3h_ líh 2 a 2I a /(h) = c0 + c, ........... pro h < a pro h > a Parametry použitého teoretického semivariogramu jsou: c0 = 2,5 ci = 7,5 a = 10,0 (dosah) Data v pěti měřených bodech mají následující souřadnice i X Y z 1 2 2 3 2 3 7 4 3 9 9 2 4 6 5 4 5 5 3 6 Pokud budeme dále značit: 47 A - matice semivariancí mezi všemi dvojicemi bodů b - vektor semivariancí mezi všemi body a bodem predikovaným X - vektor vah jednotlivých bodů 0 - tzv. Lagrangeův člen potom základní vztah pro odhad metodou krigování lze psát jako: A-Á = b Pro vlastní řešení je nutné vypočítat váhy 1, které musí splňovat podmínku ^ A = 1 Uvedený základní vztah lze vyjádřit jako soustavu rovnic: r n ľl2 ■ ■ ■ ľln 1 ľ 21 Tu ■ ■ ■ Vln 1 * = 7,1 ľni ■ 1 K 1 1 . . . 1 0 1 V tomto zápisu poslední řádek a poslední sloupec v první matici a hodnota Lagrangeova členu 0 jsou použity pro zajištění podmínky sumy vah ~^A = 1. Hodnota Lagrangeova multiplikátoru 0 také slouží pro výpočet rozptylu odhadnuté hodnoty. Uvedená soustava rovnic nám poskytne hodnoty všech vah 1 a hodnotu 0. V maticovém zápisu lze tedy psát: X A-l-b Aby bylo možné vyčíslit hodnoty semivariancí, je v prvním kroku zapotřebí vytvořit matici vzdáleností mezi datovými body: i 1 2 3 4 5 1 0, 000 5, 099 9, 899 5, 000 3, 162 2 5, 099 0, 000 6, 325 3, 606 4, 472 3 9, 899 6, 325 0, 000 5, 000 7, 211 4 5, 000 3, 606 5, 000 0, 000 2, 236 5 3, 162 4, 472 7, 211 2, 236 0, 000 Vektor vzdáleností mezi měřenými body a bodem predikovaným: i 0 1 4, 234 48 2 2, 828 3 5, 657 4 1, 000 5 2, 000 Těchto vzdáleností využijeme k výpočtu semivariancí pro sférický model semivariogramu s výše uvedenými parametry c0 , Cj, a - tedy k sestavení matice A a vektoru b: Matice A: i 1 2 3 4 5 6 1 2, 500 7, 739 9, 999 7, 656 5, 939 1, 000 2 7, 739 2, 500 8, 667 6, 381 7, 196 1, 000 3 9, 999 8, 667 2, 500 7, 656 9, 206 1, 000 4 7, 656 6, 381 7, 656 2, 500 4, 936 1, 000 5 5, 939 7, 196 9, 206 4, 936 2, 500 1, 000 6 1, 000 1, 000 1, 000 1, 000 1, 000 0, 000 Ve výše uvedené matici má řádek navíc (i=6) zajistit podmínku, že váhy budou mít sumu rovnu jedné. Vektor b: i 0 1 7, 151 2 5, 597 3 8, 815 4 3, 621 5 4, 720 6 1, 000 Inverzní matce A1: 49 1 -, 172 , 050 , 022 -, 026 , 126 , 273 2 , 050 -, 167 , 032 , 077 , 007 , 207 3 , 022 , 032 -, 111 , 066 -, 010 , 357 4 -, 026 , 077 , 066 -, 307 , 190 , 030 5 , 126 , 007 -, 010 , 190 -, 313 , 134 6 , 273 , 207 , 357 , 003 , 134 6, 873 Řešením výše uvedené soustavy rovnic lze pro jednotlivé vzdálenosti získat hodnoty vah 1: i A 1 0,0175 2 0,2281 3 0,0891 4 0,6437 5 0, 1998 6 0,1182 vypočtené hodnoty vah vypočtená hodnota Pro váhy i=l,... 5 platí, že jejich suma se rovná jedné, v posledním řádkuje hodnota Lagrangeova členu 0. Vzdálenosti měřených bodů od bodu predikovaného, již přísluší výše určené váhy: i 0 1 4, 234 2 2, 828 3 5, 657 4 1, 000 5 2, 000 Potom odhad hodnoty Z v bodě (i=0) o souřadnicích (x=5, y=5): Z(xi=0) = 0,0175*3+0,2281*4-0,0891*2+0,6437*4+0,1998*6 = Z(xi=0)=4,560 Rozptyl odhadu: 50 ae2= [0,0175*7,151+0,2281*5,597-0,0891*8,815+0,6437*3(527+0,7995*< 720]+ cp = ae2= 3,890+ 0,1182 = ae2 = 4,008 4.1 Typy krigování Na rozdíl od deterministických metod interpolace nabízí metody krigingu vedle odhadů vlastní interpolované hodnoty také odhady pravděpodobnosti výskytu těchto hodnot a dále odhady chyb predikce. Pro charakterizování jednotlivých typů krigingu předpokládejme jednoduchý model: Z(xI) = //(xI) + ^(xI) kde Z(xi) je proměnná v bodě xÍ3 která se skládá z deterministické hodnoty trendu ju(xí) a autokorelováné náhodné proměnné s(x\). Protože ve většině případů hodnotu trendu často pouze odhadujeme a neznáme ji přesně, je určitou chybou zatížena též náhodná složka. Pro tuto platí, že její průměrná chyba je rovna nule a autokorelace hodnot s(x\) efxi+h) nezávisí na aktuální pozici, ale pouze na hodnotě vzdálenosti h. Hodnota trendu může nabývat konstantní hodnoty (ji(x])= ju. V závislosti na tom, zda hodnota ju představuje konstantu či zda data obsahují trendovou složku a v závislosti na tom, za hodnota// představuje zámou hodnotu či zda ji odhadujeme definujeme různé metody (modely) krigování. Mezi základní metody krigování, kterými se provádí odhad na základě přímo naměřených hodnot patří především: • základní (ordinary) krigování s bodovým odhadem • základní (ordinary) krigování s blokovým odhadem • jednoduché (simple) krigování • univerzální krigování • pravděpodobnostní krigování • co-kriging • lognormální krigování. 4.1.1 Základní krigování (ordinary kriging) Obecný model základního krigování: Z(Xj) = //(Xj) + £(Xj) kde //je neznámá hodnota trendu. 51 \ • * • * • ......- •.............. • * ><■> • • t-1-1-1-ä-i-r D 5 10 15 20 25 30 X-Cocrflinäe Obr. 4.2 Princip základního krigování 4.1.2 Jednoduché krigování (Simple kriging) Nejjednodušší variantou krigování je tzv. jednoduché krigování (simple kriging). K výpočtu je potřebná znalost průměrné hodnoty veličiny v poli (ju) a dále předpoklad štacionarity. Obecný model jednoduchého krigování má opět tvar: Z(x1) = Ju(x1) + s(x1) kde //je známá konstanta. \ • • * * • • • Z(») * i-1-1-1-r 0 5 10 1 5 20 25 30 X-Coardincte Obr. 4.3 Princip základního krigování V tomto modelu, protože známe hodnotu /j, potom v bodech měření známe přesně také s(x). V případě základního krigování j sou obě hodnoty odhadovány. Základní korigování tedy nabízí přesnější odhad autokorelace, předpoklad znalosti //je však často nereálný. Často se však používá některá z běžných trendových funkcí, kterou se nejprve vyjádří rezidua. Trend reziduí se potom považuje za nulový a aplikuje se základní krigování. 4.1.3 Univerzální krigování (Universal kriging) Obecný model univerzálního krigování má opět tvar: Z(xi) = ju(xi) + s(xi) 52 kde ju(x) je deterministická funkce - např. polynom druhého stupně jako na přiloženém obrázku. Pokud odečteme hodnotu polynomu od originálních dat, dostaneme chybovou složku e(x), která má charakter náhodné proměnné s průměrem rovným nule. Autokorelace je modelována právě z této náhodné proměnné. Jak je patrné z obrázku, univerzální kriging je v tomto případě analogií regresní závislosti. Na rozdíl od regrese, kde složku e(x) považujeme za nekorelovanou, v případě krigování ji modelujeme jako složku autokorelovanou. Stejně jako v případě základního krigingu, vhodná dekompozice na obě výše uvedené složky ze samotných dat nelze provést. Univerzální kriging používá k popisu autokorel ováné náhodné složky semivariogramu nebo kovariance. Obr. 4.4 Princip univerzálního krigování 4.1.4 Indikátorové krigování V řadě úloh nás nezajímá, zdaje v daném místě nebo na dané ploše nej pravděpodobnější průměrná koncentrace 0,016 nebo 0,015, ale odhad pravděpodobnosti, s jakou je překročena limitní hodnota např. 0,012. Tyto úlohy řeší tzv. indikátorové krigování, které náleží do skupiny neparametrických geostatistických metod. Dále sem patří také soft kriging a pravděpodobnostní krigování. Pro tyto metody krigingu se zavádí pojem prahové hodnoty, kterou lze spojitou proměnnou převést na binární (viz. obr.) Indikátorové krigování předpokládá model ve tvaru: /(*,.) = // + *>(*,.) kde ju neznámá konstanta, e(x) autokorel ováná náhodná proměnná a I(x) je buďto přímo zjištěná binární proměnná a nebo binární proměnná, kterou obdržíme ze spojitých dat prahováním. Jinak se model neliší od základního krigování. Protože indikátorová proměnná nabývá hodnot 0 nebo 1, potom interpolované hodnoty budou nabývat hodnot od 0 do 1 a lze je interpretovat jako pravděpodobnost, že proměnná nabývá hodnoty 1. Pokud bylo pro vytvoření indikátorové proměnné použito postupu prahování, potom lze vytvořenou mapu interpretovat jako pravděpodobnost překročení prahové hodnoty. 53 > •ä Values abtľve iteeshold become l's. ♦ » Threshold * Valuta below threshold become O's r •Ji 25 —r- Kl 10 15 X-C □ordinate Obr. 4.5 Princip prahování hodnot proměnné \ □□□□□□□□□□□□□□□□□□ □ □□□□□□□ ..................i.................................r............................. o s 10 15 20 2fl 30 X-Coordinate 1 > & IT: Obr. 4.6 Princip indikátorového krigování Uvedený princip lze obecně rozšířit a použitím dvou více hodnot vytvořit např. dvě indikátorové proměnné (viz. dále - co-kriging). 4.1.5 Pravděpodobnostní krigování (probalility kriging) Pravděpodobnostní krigování předpokládá model ve tvaru: l(x) = I(Z(x) >cl) = jul+sl (x) Z(x) = Ju2+s2(x) kde jui a p.2 jsou neznámé konstanty. I(x) je binární proměnná vytvořená indikátorovým prahováním (I(Z(x) > ci). V tomto případě dostáváme dvě náhodné chyby &i(x) a &2(x). Cíle pravděpodobnostního krigování j sou stejné jako u krigování indikátorového, j sou však dosaženy využitím konceptu co-krigingu. Na obrázku 7 má datový bod Z(u=9) hodnotu indikátorové proměnné l(u)=0 a bod Z(x=10) hodnotu I(x)=l. Pokud bychom chtěli predikovat hodnotu v polovině vzdálenosti mezi oběma body - na x-ové souřadnici 9,5, potom použitím modelu indikátorového krigování bychom obdrželi hodnotu 0,5. Z obrázku je však patrné, že datový bod Z(x) je nepatrně nad hodnotou 54 prahu, naopak bod Z(u) je výrazně pod prahovou hodnotou. Je tedy reálné předpokládat, že predikovaná proměnná v bodě 9,5 bude méně než 0,5. > OJ ôi OJ □□□□□□□□□□□□□□□□□n • • • . * • • * • • • • • • XZ(x) ^Z(u) +—1(11) Ct • 1.0 0-8 0.6 0.4 0.7 0 0 10 15 20 X-C o ordinate 30 Obr. 4.7 Princip pravděpodobnostního krigování Pravděpodobnostní krigování se tedy snaží využít vedle indikátorové proměnné ještě další extra informace v původních datech. Nevýhodou pravděpodobnostního krigování je nutnost provádět odhady jako autokorelací pro jednotlivé proměnné, tak křížových korelací mezi mini. Dalšími odhady neznámých autokorelací se vnáší do výsledného modelu větší míra nejistoty. 4.1.6 Nelineární kriging (log-normal) Pokud nemají vstupní data normální rozdělení, je nutné je před vlastní interpolací transformovat. Nejběžnější je transformace lognormální. Originální data jsou transformována na přirozený logaritmus o základu 10. Tedy modelování variogramu a interpolace probíhá s proměnou y(u): y(u) = ln z(u) Predikované hodnoty je poté nutno transformovat nazpět, což může působit problémy (viz. Borrough et. al. 1992) a jako alternativa se nabízí indikátorový kriging Pro některá FG data, která vykazují rozdělení s kladnou asymetrií, je však lognormální transformace výhodná (např. obsah chemických látek v půdě). 4.1.7 Kriging s využitím externí informace K interpolaci kromě hodnot vlastní interpolované proměnné lze využít například: 1. vhodnou stratifikaci dat (stratifikovaný kriging) 2. hodnoty jiné proměnné, která koreluje s původní a kterou lze snadno měřit ve větším počtu bodů (např. výškové poměry) - co kriging 3. fyzikální či empiricky sestavený model, který podmiňuje rozložení hodnot studované proměnné 55 Stratifikovaný kriging spočívá v rozdělení oblasti na subregiony. Předpokládá dostatečný počet bodů pro výpočet hodnot variogramu. Může dávat vhodnější odhady, je však nutné řešit oblasti na styku subregionů. Např. obsah znečišťujících látek podle oblastí zaplavovaných podél vodního roku s různou frekvencí. 4.1.8 Co-kríging Máme dvě proměnné zi a Z2, které vykazují prostorovou korelaci. Pak lze využít hodnot proměnné z2, k interpolaci hodnot proměnné z\. Tento koncept je vhodný zvláště v případech, kdy je proměnná z2 snáze získatelná a rozšiřitelný i na více než dvě proměnné. Přitom pro přesnější odhady se používá jak autokorelace jednotlivých proměnných, tak vzájemné (cross) korelace všech použitých proměnných. Základní co-kriging využívá následujících modelů: kde jui a//2 jsou neznámé konstanty. Dále dostáváme dvě náhodné chyby ei(x) a e2(x). Základní co-kriging odhaduje hodnotu proměnné Zi(xo) stejně jako základní krigování, ovšem navíc využívá kovariance s hodnotu Z2(x). 0 5 10 1S 20 25 30 X-Coordinate Obr. 4.8 Princip co-krigingu Z obrázku je patrné, že data Zi a Z2 se jeví jako nekorelovaná. Dále pokud Zi je pod průměrem fi\ , potom Z2 je často nad průměrem //2 a naopak. Tedy Zi a Z2 vykazují negativní cross korelaci. Vedle základního co-krigingu jsou dalšími variantami např. jednoduchý, univerzální, indikátorový či pravděpodobnostní co-kriging. 4.1.9 Blokový odhad při základním krigování (Block kriging) Lokální (bodový) odhad metodou krigingu lze určitým způsobem vztáhnout k ploše či objemu v prostoru interpolovaných dat. Mnoho přírodních jevů vykazuje značnou variabilitu a výsledkem bodového odhadu může být mapa obsahující značný počet ostrých vrcholů a depresí. Tento efekt lze potlačit tak, že modifikujeme výše uvedené rovnice a odhadneme 56 průměrnou hodnotu z(B) proměnné z pro jistou plochu či objem B (viz. obr). Tato modifikace je vhodná, pokud výsledkem interpolacemi být struktura pravidelných buněk (grid). Z4 Obr. 4.9 Princip blokového krigování Průměrná hodnota z pro blok B z(B) = f Z(x)dX J plocha _ B bude odhadnuta z výrazu: i=l Kde stejně jako u bodového odhaduje suma všech vah Xi rovna jedné. Minimální rozptyl nyní bude: o-2(B) = |]AIř(xI,B) + ^-ř(B,B) i=l a získáme ho, když n 2]\ľ{\, Xj) + = ř(Xj, B) pro všechna j. i=l Rozptyly odhadů pro blokový kriging jsou daleko menší než pro bodový kriging. Výsledný interpolovaný povrch je obecně více shlazený a neobsahuje takové množství lokálních extrémů. Blokové korigování je aproximující metodou. 4.2 Hodnocení a verifikace modelů Krigování jako interpolační metoda umožňuje pro každý interpolovaný bod odhadnout potenciální velikost chyby odhadu. Vedle map predikovaných hodnot tak lze především konstruovat mapy hodnot j£i_=<----• 1 11 ................................i.....• ^J»'.i»r* *................. a 0 89 • ..... » ■ i..... 0.68 rV^ A n 4b' I 046 0 68 0 89 1 10 1.31 1.52 1 74 Measured, 101 Regression function 0 778" x + 0.022 Obr. 4.12 Korelační pole měřených a predikovaných hodnot Chybový graf (Error plot) - stejný jako předchozí, jsou však vynášeny hodnoty rozdílů mezi měřenými a predikovanými hodnotami Standardizovaný chybový graf (Standardized Error) - hodnoty rozdílů mezi měřenými a predikovanými hodnotami j sou děleny odhadnutou směrodatnou chybou krigování. V případě nulové autokorelace budou všechny predikované hodnoty stejné - budou odpovídat průměru a proložená přímka bude mít horizontální průběh. V případě prostorové autokorelace a vhodného modelu krigingu bude proložená přímka totožná s diagonálou a navíc body korelačního pole budou vykazovat malé odchylky od diagonálního směru. Q-Q graf-znázorňuje graf kvantilů rozdílů mezi měřenými a predikovanými hodnotami dělenými odhadnutou směrodatnou chybou krigování a odpovídajících kvantilů normovaného normálního rozdělení. V případě, že odchylky měřených a odhadnutých hodnot mají normální rozdělení, potom se body v korelačním poli přimykají k přímce (viz. obr.) Predicted | Error | Standardized Error QQPÍot ,| S 3.06 1 209 S 112 1 015 1 -0 82 8 -179 -2.76 ----------- i .......................»..»_«.• ___r9'~~ -2. 3 -1.86 -0.93 -0.00 0.93 1.86 2.80 Normal Value Obr. 4.13 Příklad Q-Q grafu 61 4.3 Interpretace statistických charakteristik k hodnocení vhodnosti modelu: • Požadavek nestrannosti odhadu - unbiased - průměrná chyba odhadu a standardizovaná průměrná chyba odhadu by se měly blížit k nule: o MPE -» 0 o MSPE^O • Požadavek minimálních chyb - aby predikované hodnoty byly co nej blíže hodnotám měřeným. Čím menší bude hodnota RMSPE, tím lepší model - tedy tuto podmínku lze použít k porovnání vhodnosti více modelů. o RMSPE -» min. • Požadavek vhodné variability předikovaných dat - variabilita předikovaných hodnot je určována z hodnot měřených. Je tedy důležité, aby i variabilita interpolací vypočtených hodnot byla vhodná: o ASE « RMSPE - vhodný model (vhodná variabilita předikovaných hodnot) o ASE > RMSPE - máš model nadhodnocuje variabilitu odhadnutých hodnot o ASE < RMSPE - máš model podhodnocuje variabilitu odhadnutých hodnot V případě značného podílu šumové složky (např. v důsledku chyb v měření) či v případě značně komplexního povrchu nedává kriging lepší výsledky než jiné interpolátory. Na rozdíl o jiných metod kriging nabízí objektivní, a priori metodu odhadu vhodného okolí pro vlastní interpolaci. Řeší tedy otázku počtu bodů v okolí daného bodu, otázku velikosti a tvaru tohoto okolí. V případě existence bariér (náhlých skoků v hodnotách interpolovaného povrchu nedává kriging dobré výsledky a je nutné jej rozdělit na elementární části neobsahující bariéry. 5. Modelování prostorového uspořádání bodů Deskripce bodů pomocí měr úrovně a variability je jen prvním krokem analýzy. V případě prostorové analýzy nás v druhém kroku zajímají body s ohledem na jejich prostorové rozmístění (strukturu - pattern). Rozmístění boduje výsledkem určitých procesů a podmínek - např. lokace měst je výsledkem působení faktorů jako je reliéf, přírodní zdroje, komunikace, obdobně výskyt rostlinných druhů, atd. Cílem studia prostorového rozmístění je zjistit, jak daleko má konkrétní rozmístění objektů k rozmístění teoretickému, (např. teorie centrálních míst - teoretický vzorec - šestiúhelníky). To nám umožňuje jednak porovnávat rozmístění objektů pro různé prostorové jednotky (kategorie landuse, půdní typy, okresy, státy, atd.), jednak studovat dynamiku změn v rámci jedné jednotky (studium dynamiky). Statisticky prokázaný výskyt určitého prostorového uspořádání (shlukového či pravidelného vzorku) může být základem pro zjišťování příčin, které vedly k pozorovanému uspořádání 62 5.1.1 Statistická deskripce prostorových vzorů bodových prvků ******* ******* ****** li ******* * ****** ******* ******* □nnn □□□c □ne □□□ □□□□ □□□ mr □□□c □□□□ IIIII i n II II Obr. 5.1 Základní typy prostorového uspořádání bodů (1. sloupec), linií (2. sloupec) a ploch (3. sloupec).Typy uspořádání: shlukové (1. řádek), pravidelné (2. řádek), náhodné (3. řádek) Rozlišujeme tři základní typy prostorového uspořádání bodů: • Shlukové (Clustered) • Pravidelné (Regular) • Náhodné (Random) 5.1.2 Základní metody statistického popisu prostorového uspořádání bodů: • Analýza kvadrátů - testujeme, zda rozmístění bodů v ploše je náhodné či nikoliv. • Metoda nejbližšího souseda - porovnává průměrnou vzdálenost mezi nejbližšími sousedy pole bodů vzhledem k teoretickému rozmístění. • Metody prostorové autokorelace - měří, jak podobné či nepodobné jsou hodnoty atributů sousedních bodů. 5.1.3 Problém měřítka, rozsahu studované oblasti a kartografické projekce Měřítko - je nutné vhodně zvolit tak, aby studovaný jev mohl být prezentován body v prostoru. Rozsah studované oblasti - v závislosti na zvolené oblasti (často vymezené administrativními hranicemi) se mění jak vzdálenosti mezi jednotlivými body, tak také charakteristiky jejich prostorového uspořádání (Obr. 5.2). 63 Obr. 5.2 Vliv velikosti studované oblasti na prostorové uspořádání bodů Kartografickou projekci je nutno vhodně zvolit podle účelu (viz. Analýza kvadrátů). Projekcí se mění tvar, vzdálenosti, vzájemná poloha objektů (viz. Obr. 5.3). Čím větší studovaná oblast, tím větší bude role zvolené projekce. Obr. 5.3 Vliv kartografické projekce na tvar studované oblasti 5.2 Analýza kvadrátů (QUADRAT ANALYSIS) Metoda pro detekci prostorového uspořádání bodů. Je založena na hodnocení změn hustoty bodů v prostoru. Je porovnáváno, zda rozmístění bodů v prostoru je náhodné, či má blíže k uspořádání shlukovému či pravidelnému. _|n|* VÍ Ottcourtty.shp CD -í -> r Obr. 5.4 Analýza kvadrátů - pravidelné rozmístění buněk 64 Postup analýzy spočívá v rozdělení studované plochy pravidelnou sítí na buňky a je zjištěn počet bodů v každé buňce. Následně je analyzováno rozdělení četností buněk s určitým počtem bodů. Toto rozdělení je porovnáváno s náhodným rozdělením četností. Buňky se označují jako kvadráty a nemusí jít o čtverce, ale např. i o kruhy či šestiúhelníky. Tvar buněk většinou vychází z empirie. V rámci jedné analýzy však tvar a velikost buněk musí být konstantní. Extrémně shlukové uspořádání - většina bodů v jedné či několika málo buňkách Extrémně pravidelné - ve všech buňkách přibližně stejně Uvedenou metodu lze využít také tak, že se buňky stejné velikosti náhodně rozmístí po studované ploše. Obr. 5.5 Analýza kvadrátů - náhodné rozmístění buněk Citlivou stránkou metody je volba velikosti kvadrátů. Optimální velikost kvadrátů (QS) lze získat z následujícího vztahu: QS=^ n kde A je plocha studované oblasti a n počet analyzovaných bodů. Velikost strany vhodného kvadrátu j e potom: Získané rozložení četností bodů v kvadrátech (empirické) je porovnáváno s náhodným rozložením (teoretickým). Vhodným testem je např. K-S test. Testem můžeme kvantifikovat rozdíl empirického a teoretického (shlukové, pravidelné, náhodné) rozdělení bodů v ploše. 5.2.1 Praktický postup testování výsledků analýzy kvadrátů: Formulujeme nulovou hypotézu - neexistuje statistiky významný rozdíl Qe-\i rozdíl malý, může být výsledkem náhody, čím je větší, s tím větší pravděpodobností náhodný není, aleje statistiky významný). Zvolíme hladinu významnosti a = 0,05 65 Vypočítáme kumulované četnosti Vypočteme testovací kritérium: Vypočteme kritickou hodnotu D = max Q - B D =í# V m kde m je počet kvadrátů. V případě porovnávání dvou výběrů o různém počtu členů ml a m2 se kritická hodnota vypočte následovně: D„=u6Jä±ä: \ m, -rrij Je-li vypočtená hodnota D větší než kritická hodnota Da, potom rozdíl mezi oběma uspořádáními je statisticky významný. Pozorované rozložení bodů můžeme také porovnávat s rozložením náhodně generovaným (např. podle určitého teoretického rozdělení). Často se využívá rozdělení Poissonovo (Poisson random process) Poissonovo rozdělení je určeno především průměrnou frekvencí výskytu (k) v jednotlivých jednotkách (kvadrátech), kde A = ^/^ při m kvadrátech a n bodech v prostoru. Je-li x počet bodů v kvadrátu, potom pravděpodobnost výskytu x bodů v kvadrátu podle Poissonova rozdělení je definována: e a P(x) = —- x! Z uvedeného vztahu můžeme pro různá x vypočítat pravděpodobnost rozložení bodů, které budou mít Poissonovo (náhodné) rozdělení. Hodnoty pravděpodobnosti lze zjistit i zkráceným výpočtem. Je-li x=0, potom p(0) = e~Ä a pravděpodobnosti pro následná x můžeme určit z p(0), obecně: Ä P(x)= p(x-l)* — x Je-li x= 1, potom p(x-l) = p(0) atd. Vedle K-S testu můžeme k hodnocení rozdělení bodů v kvadrátech použít také vlastností Poissonova rozdělení - především hodnoty průměru a rozptylu Poissonova rozdělení, pro které platí, že se rovnají hodnotě QC). Jinými slovy bude-li distribuce bodů v prostoru 66 generována náhodným procesem, potom toto rozdělení má stejný průměr a rozptyl. Tedy jejich poměr se bude blížit jedné. Postup: Vypočteme hodnoty průměru a rozptylu pro četnosti bodů v kvadrátech a hodnoty dáme do poměru. Hodnotu porovnáme s 1. Rozdíl lze dále standardizovat (vyjádřit v násobcích směrodatné odchylky). Vyjde-li hodnota větší než 1,96, potom je rozdíl statisticky významný na hladině a = 0,05. Test založený na poměru průměru a rozptyluje silnější než K-S test, lze ho však použít pouze v případě, že předpokládáme Poissonovo rozdělení studované množiny bodů. Pozorované rozdělení bodů lze porovnávat i vůči jiným teoretickým rozdělením (např. negativní gamma či negativní binomické). Omezení analýzy kvadrátů: Obr. 5.6 Analýza kvadrátů neřeší otázku rozložení bodů uvnitř kvadrátů 5.3 Analýza nejbližšího souseda (NEAREST NEIGHBOUR ANALYSIS) Metoda analýzy kvadrátů je založena na konceptu hustoty (počet bodů v ploše). Metoda analýzy nejbližšího souseda je naopak založena na konceptu vzdálenosti (spacing - plocha připadající na bod). Metoda analýzy nejbližšího souseda je založena na porovnání pozorované průměrné vzdálenosti mezi nejbližšími sousedy a této průměrné vzdálenosti u známého vzorku (pattern). Pozorovaná průměrná vzdálenost mezi nejbližšími sousedy může být větší či menší než vzdálenost při náhodném rozmístění bodů. Obr. 5.7 Analýza nejbližšího souseda - pravidelné uspořádání bodů 67 Homogenní oblast - nejvíce uniformní vzorek - body v ploše tvoří středy pravidelných šestiúhelníků. Body tvoří trojúhelníkovou mřížku. Za této konfigurace bude vzdálenost mezi body rovna výrazu kde A je plocha a n počet bodů v ploše. V reálné situaci tvoří geografické rozložení bodů výjimečně pravidelný vzorek. K testování, zda má určité rozložení bodů v ploše jistý vzorek lze využít R statistiku (R - randomness). Určí se jako poměr mezi pozorovanou a očekávanou průměrnou vzdáleností nejbližších sousedů v určité oblasti: Hodnotu r0bs zjistíme tak, že určíme vzdálenost mezi daným bodem a všemi jeho sousedy. Dále najdeme nejkratší vzdálenost - tedy nejbližšího souseda. Tento proces se opakuje pro všechny body. Ze všech nej kratších vzdáleností se vypočte průměr. Pro teoretické - náhodné - rozložení se průměrná vzdálenost nejbližšího souseda vypočte podle vzorce: 1 exp Čím je hodnota R < 1, tím více se prostorové rozložení bodů blíží rozložení shlukovému (robs"^ rexp). Čím je hodnota R > 1, tím více se prostorové rozložení bodů blíží rozložení pravidelnému (ľobs ľexp). R=0 R=0.51 R=1.0 R=1.48 R=1.90 SHLUKOVÉ PRAVIDELNE Obr. 5.8 Škála hodnot R statistiky • R = 2,149 • R = 0 • R= 1 zcela shlukové náhodné zcela pravidelné 68 Je-li R=0, vzdálenosti jsou 0, všechny body mají stejnou polohu. Jinou z možností, jak porovnat rozdíl mezi pozorovanou a očekávanou vzdáleností nejbližšího souseda je porovnat tuto diferenci s tzv. směrodatnou chybou (Standard Error - SEr) Směrodatná chyba popisuje pravděpodobnost, že jakýkoliv rozdíl dvou hodnot je výsledkem náhodných vlivů. Je-li tedy zjištěná diference malá ve srovnání s SE, potom rozdíl není statisticky významný a naopak. Použití směrodatné chyby SE vychází z vlastností normálního rozdělení, pro které platí následující: Je-li mezi pozorovanými populacemi rozdíl a jeho velikost náleží do intervalu (-lSEr; + lSEr), potom existuje 68 % šance, že tento rozdíl je náhodný - tedy nevýznamný: Pravděpodobnost (<68%) = f-lSE,; +lSEr) Za statisticky významný považujeme rozdíl, který můžeme obdržet v 5 případech ze sta - tedy s pravděpodobností 5 %, a=0,05. Vyjádřeno v násobcích směrodatné chyby - rozdíl mezi dvěma populacemi povařujeme za statisticky významný, jestliže je menší než -l,96SEr a nebo větší než + l,96SEr: Pravděpodobnost (<95%) = f-l,96SEr; + l,96SEr) Výpočet směrodatné chyby pro pozorované vzdálenosti bodů: cc 0,26136 ODr — -, VnVl Pomocí směrodatné chyby lze vypočítat standardizovanou hodnotu (Z-score): r — r ry _ obs exp R~ SEr Je-li tedy Zr< -1,96 či Zr> 1,96 potom vypočtený rozdíl mezi pozorovaným a náhodným uspořádáním je statisticky významný - tedy není náhodný a naopak. Nelze spoléhat na vizuální srovnání prostorového rozložení ani na vypočtenou hodnotu R. Ta by měla být doplněna hodnotou ZR pro ověření statistické významnosti pozorovaného rozdílu. Metoda analýzy nejbližšího souseda může být rozšířena na analýzu nejbližších sousedů druhého, třetího a vyšších řádů. Například u obr. 2.6 dokumentujícího nevýhody kvadrantové analýzy by až analýza nejbližšího souseda druhého řádu odhalila, že se obě uspořádání výrazně liší. Na obrázku vlevo je R-statistika druhého řádu velká, na obrázku vpravo naopak malá. Použití analýzy nejbližšího souseda rozdílných řádů může odhalit heterogenity v uspořádání bodů na rozdílných prostorových úrovních. 69 Problémy spojené s metodou analýzy nejbližšího souseda: výsledky jsou vysoce citlivé k měřítku (lokální vs. regionální) a vymezení zpracovávané oblasti. V závislosti na studovaném jevu by měla být věnována pozornost také vymezení studované plochy (administrativní či přirozené hranice). 5.4 Prostorová autokorelace (SPATIAL AUTOCORRELATION) Jak analýza kvadrátů, tak analýza vzdálenosti nejbližšího souseda pracují pouze s polohou bodů. Nerozlišují body podle hodnot jejich atributů. Oba parametry (polohu i atributy) hodnotí prostorová autokorelace - je tedy metodou vhodnější. Východiska prostorové autokorelace: Většina jevů se v prostoru mění spojitě. Blízké body budou mít i podobné hodnoty studovaného jevu a naopak. (First law of geography - Tobler, 1970) Koeficient prostorové autokorelace - uvažuje polohu bodů (vzájemnou vzdálenost) a hodnotí rozdílnost hodnot atributů bodů v prostoru. Mezi nej používanější koeficienty prostorové autokorelace náleží Gearyho poměr C (Geary's Ratio) a Moranův index I (Moran's I). Lze jich využít pro intervalová a poměrová data. Dále používaná notace: Cjj - podobnost atributu v bodě i a j • Wjj - vzdálenost bodu i a j. wu = 0 pro všechny body • xi - hodnota studovaného atributu v bodě i • n - počet bodů ve vyšetřovaném vzorku Obě míry prostorové autokorelace kombinují v jednom výrazu míry podobnosti atributů i míry podobnosti polohy - tento výraz je potom východiskem pro definování dalších vztahů: Koeficient prostorové autokorelace SAC (spatial autocorrelation coefficient) je úměrný vážené míře podobnosti atributů bodů - obecně: n n IIv^ SAC^^lii- II i=i j=i W" y V případě Gearyho poměru se podobnost hodnot atributu mezi dvěma body vypočte podle následujícího vztahu: qj=(jq-xj)2 70 Gearyho poměr C se tedy vyjádří jako: n n n n c i=l j=l_ _ i=l j=l i=i j=i i=i j=i kde g2 je rozptyl hodnot atributu x s průměrem x Z(^-x)2 .2 _ _i = l_ (n-1) V případě hodnoty Moranova indexu I se podobnost hodnot atributu v bodech i a j vyjádří následovně: Cg = (Xj - x) • (Xj - x) Moranův index I je potom určen: n n n n ZZV^j SXwij • (Xi - x) • (Xj - x) j _ i=l J=l _ i = l J=l n n n n í=i j=i i=i j=i kde s je v tomto případě výběrový rozptyl: ,2 _ i=l Ve výše uvedených vzorcích lze všechny neznámé přímo určit z hodnot atributů bodů. Jedinou doposud nedefinovanou neznámou zůstává míra podobnosti (blízkosti) polohy bodů i a j, tedy hodnota Wy. Ta se běžně uvažuje jako inverzní hodnota vzdálenosti těchto bodů. Tedy podle výše uvedených předpokladů dáváme malou váhu hodně vzdáleným bodům a velkou váhu hodně vzdáleným bodům, tedy: 5.4.1.1 Wj = V, Rozdíly mezi oběma indexy jsou dány způsobem výpočtu rozdílů mezi hodnotami atributu. Obor hodnot, kterých mohu oba indexy nabývat se tedy také liší, jak uvádí následující tabulka: Prostorové uspořádání Gearyho poměr C Moranuv index I 71 Shlukové uspořádání, sousední body vykazují podobné hodnoty 0E(l) Náhodné uspořádání, body nevykazují znaky podobnosti C~ 1 I = E(l) Pravidelné uspořádání, sousední body vykazují rozdílné charakteristiky 1 < C < 2 I < E(l) kde E(I) = (-l)/(n-l) 5.4.2 Předpoklad náhodnosti a předpoklad normality Při studiu prostorového uspořádání, můžeme předpokládat dva základní způsoby, kterými jsou atributy přiřazeny j ednotlivým bodům. Předpoklad náhodnosti (randomization, nonfree sampling) - předpokládáme, že hodnoty atributů v bodech představují pouze jednu z možných variant uspořádání při použití stejné množiny hodnot. Alternativně můžeme předpokládat, že hodnoty atributů v množině studovaných bodů jsou pouze jednou z nekonečného množství možností. Každá hodnota je nezávislá na hodnotách jiných v množině bodů - předpoklad normality (normality, free sampling). Předpoklad normality dovoluje nahrazení hodnot pozorování na rozdíl od předpokladu náhodnosti. 5.4.3 Určení odhadů očekávaných hodnot Výše uvedené předpoklady náhodnosti (R) a normality (N) ovlivňují způsob výpočtu očekávaných (e - expected) hodnot i hodnot rozptylu. Očekávané hodnoty indexů a hodnoty rozptylů potřebujeme pro testování, zda se vypočtené hodnoty indexů C a I statisticky významně liší od náhodného uspořádání. Odhad očekávaných hodnot pro náhodné uspořádání (random pattern) a rozptyly pro Gearyho poměr C: VARN(C) EN(C) = 1 ER(C) = 1 [(2S1 + S2)(n-1)-4W2] 2(n + l)W2 VARR(C) kde (n - X)SX [n2 - 3n + 3 - (n - l)k] (n - 1)S2 [n2 + 3n - 6 - (n2 - n + 2)k] W2[n2 -3 -jn-lfk] n{n - 2)(n - 3)WZ An{n-T)(n-V)W n{n - 2)(n - 3)WZ w = ZZ i=l j=l 72 ZV (w. +w..)2 s2 = Z(wí. + wí)2 i=l T" (x,-x)4 ín V Z(^-x)2 U=i J Očekávané hodnoty Moranova indexu I a hodnoty rozptylu se pro náhodné uspořádání vypočtou obdobně: EN(l)=ER(l) = — n-1 VARN(I)- . . N W2(n2-1) VAR (J)_n\{n2-3n + 3)si-nS2+3W2] k\{n2-nS2 +3W2] r ■» rK ' (n-\)(n-2)(n-3)W2 (n-\)(n-2)(n-3)W2 L R J Máme-li vypočteny očekávané hodnoty indexů a jejich rozptyly, můžeme vyjádřit standardizované hodnoty (Z-skore) Z = I-E(I) VAR(I) nebo C-E(C) Z VAR(C) Pro hodnoty Z pak mohou být použity stejné kritické hodnoty, tedy na hladině významnosti a=0,05: -1,96 < Z< +1,96 73 Obr. 5.9 Příklad výpočtu měr prostorové autokorelace Interpretace hodnot koeficientů prostorové autokorelace: Pokud zjištěné hodnoty z-skóre padnou vně intervalu (-1,96 ; +1,96), potom se prostorové uspořádání bodů statisticky významně liší (na hladině 5 %) od uspořádání náhodného. 5.4.4 Alternativy výpočtu: V uvedených vztazích lze modifikovat výrazy pro vyjádření podobnosti polohy. Například hodnoty Wy mohou nabývat binárních hodnot 0, 1 podle toho, zda jde o body sousední či nikoliv (viz. např. teorie nodálních regionů, kdejako sousední body považujeme centroidy regionů, které obklopují daný region. Modifikovat lze také váhy vzdálenosti bodů výrazem: kde koeficient b může nabývat různých hodnot v závislosti na povaze studovaného problému (vzdálenost měřená dosažitelností autem a letadlem je jiná). Hodnota b je často rovna 2. Uvedených koeficientů prostorové autokorelace lze využít pro výpočet podobnosti mezi polygony (viz. dále). 6. Statistická analýza liniových prvků Linie mohou na mapách reprezentovat dva příbuzné objekty: • Vlastní linie - reprezentují a lokalizují skutečně lineární geografické fenomény (řeky, silnice, potrubí) • Hrany - rozdělují plochy a povrchy (hraniční linie, lomové linie). Hrany nemají šířku. 74 Problémy prezentace „přirozených linií" v prostředí GIS jsou spojeny především s procesy generalizace a zjednodušení průběhu. Linie je prezentována jako spojnice posloupnosti lomových bodů, mezi lomovými body je rovná. Problém měření vzdáleností - Někdy se místo měření vzdálenosti v délkových j ednotkách používá cestovní čas a dopravní náklady. Pro analýzu linií jsou vedle délky významné také atributy jako orientace, směr či spojení. Existence spojení mezi soustavou bodů, které tvoří linii, znamená, že lokace (body) na sobě nejsou nezávislé, ale jsou spojené v určitém směru. Body spojené v určitém pořadí musí zachovávat tuto posloupnost. Obr. 6.1 Liniové prvky na digitální mapě - prosté linie, trajektorie, síť Linie mohou v GIS vystupovat na třech úrovních, které představují jistou hierarchii (Obr. 6.1): 1. „Prosté" linie - např. zlomy - lze určit jen délku a orientaci. Může existovat jako jednoduchá spojnice dvou bodů či jako „řetězec" 2. „Trajektorie" - vektor pole větru - lze určit velikost (délku), orientaci a směr 3. Sítě - dopravní sítě, říční síť - lze určit prostorové uspořádání - topologické vztahy, konektivitu, dostupnost, ... Geometrické charakteristiky - linie může být prezentována jako: • Jednoduchá spojnice - pouze dvou bodů (koncový a počáteční - délka je Euklidovská vzdálenost • Posloupnost několika liniových segmentů - řetězec Příklady analýzy prostorových vazeb liniových prvků: • analýza převládající orientace, průměrné délky spoje, • charakterizování liniových vzorků - „uspořádání sítí" • dopravní dostupnost • gravitační modely • hledání optimální trasy 75 6.1 Prostorové atributy liniových prvků Délka linie může být definována jako: • přímá vzdálenost (vypočtená z Pythagorovy věty) • „skutečná" vzdálenost (součet přímých vzdáleností jednotlivých segmentů) Orientace linie - orientace neurčuje směr (např. JV = SZ) - orientace zlomů, ulic. Nemá smysl otázka odkud - kam? Směr linie - typicky - vektor pole větru 6.1.1 Topologie (sítí) Výše uvedené atributy linií lze vyjádřit i pro jednotlivé segmenty sítě či pro celou síť jako celek (průměrná délka sítě, převládající orientace či směr segmentů sítě). Vedle toho jsou pro charakterizování sítí důležité atributy popisující jejich strukturu a uspořádání jako celek a dále popisují vztahy segmentů uvnitř sítě (topologii). Obr. 6.2 Příklad sítě Tabulka 6.1 Matice konektivity ID 1 2 3 4 5 6 7 8 9 10 0 1 1 0 0 0 1 0 0 1 2 1 0 1 0 0 0 0 0 0 0 3 1 1 0 1 1 0 0 0 0 0 4 0 0 1 0 1 1 0 0 0 0 5 0 0 1 1 0 1 0 0 0 0 6 0 0 0 1 1 0 0 0 0 0 76 7 1 0 0 0 0 0 0 1 1 1 80000001010 90000001 100 10 1 000001 000 Základním topologickým aspektem sítě je způsob propojení jednotlivých segmentů - tedy její konektivita. Tradičním nástrojem používaným k charakterizování konektivity je matice konektivity. Je to matice čtvercová, binární, symetrická o n řádcích (sloupcích), kde n je počet segmentů sítě. Jednička v matici značí, že dva příslušné segmenty jsou bezprostředně spojeny. Na hlavní diagonále matice jsou nuly. 6.2 Směrová statistika (Directional statistics) Topologii sítě lze charakterizovat jednoduchými mírami. Takovou je např. poměr mezi skutečnou délkou linie a spojnicí počátečního a koncového bodu. Tato charakteristika se určuje křivost linie (sinusoity). Čím větší číslo, tím větší křivost. Směr linie - vizuální hodnocení směru linií lze provést přidáním šipek. Např. u pole větru je možné odhalit strukturu proudění v celé oblasti. 6.3 Směrový průměr (directional mean). Využití klasických měr popisné statistiky pro charakterizování směru a orientace linií je nevhodné (viz. obr. 3.3). Obr. 6.3 Problém popisné statistiky při určování charakteristik směru linie Jak je patrné z obrázku, aritmetický průměr dvou vektorů s úhly 45 a 315 stupňů dává 180 (jižní směr), avšak měl by být 0 stupňů (severní směr). Průměrný směr je však nutné určit vektorovým součtem či tzv. směrovým průměrem (directional mean). Protože pracuje se směrem (úhlem) a ne s délkou, je možné ho prezentovat na základě jednotkových vektorů. Vektorovým součtem - přidáním počátku druhého vektoru na konec prvního dostaneme směrový průměr. 77 Obr. 6.4 Koncept směrového průměru Směr výsledného vektoru lze získat také z následujícího vztahu: tan#„ oy ox kde oy je suma délek vektorů ve směru osy y a ox suma délek vektorů ve směru osy x. Protože všechny vektory jsou jednotkové, délka ve směru osy y je v podstatě sin úhlu a délka na ose x. je cosinus úhlu. Potom, jsou-li vektory označeny a, b, c a odpovídající úhly 8a, 6>b, 6C, potom: sin 0 + sin 6, + sin G tan6>--1-h-c- cos 9a + cos 9h + cos 9C Obecně, máme-li n vektorů v a úhel vektoru v od osy x je 6>v, výsledný vektor OR má úhel 6>R, měřený proti směru hodinových ručiček od osy x: tan#„ 2>inflv 2>os6>v což je tedy tangenta úhlu výsledného vektoru. Směrový průměr je potom ar etan z výše uvedeného výrazu. Výsledná hodnota směrového průměru musí zohledňovat specifika jednotlivých kvadrantů, jak uvádí následující pravidla: 1. čitatel i jmenovatel v tan 6>R j sou oba kladné - není nutná žádná úprava (vektor leží v 1. kvadrantu) 2. čitatel je kladný jmenovatel záporný - směrový průměr bude 180 - 6>R, (vektor leží v 2. kvadrantu) 3. čitatel i jmenovatel v tan 8R jsou oba záporné - směrový průměr bude 180 + 6>R, (vektor leží v 3. kvadrantu) 4. čitatel je záporný, jmenovatel kladný - směrový průměr bude 360 - 6>R, (vektor leží v 4. kvadrantu Praktický výpočet spočívá v určení sin a cos úhlů všech vektorů. Určí se jejich sumy a vytvoří poměr, který je tangentou výsledného úhlu. Směrový průměr je potom aretan. 78 6.4 Směrový rozptyl (Circular variance) Stejně jako v případě klasické popisné statistiky je charakterizování souboru prvků pouze měrou úrovně, kterou je výše uvedený směrový průměr, je často nedostatečné a může být i zavádějící. Např. pokud dva vektory budou svírat úhel 180 stupňů. Proto je nutné použít i měr variability (rozptylu). Pokud dáme dohromady vektory podobného směru, výsledný vektor bude relativně dlouhý. Jeho délka se bude blížit n, pokud bude n jednotkových vektorů. Naproti tomu, pokud dáme dohromady vektory opačného či značně rozdílného směru, výsledný vektor bude významně menší než n. Tedy délku výsledného vektoru můžeme použít jako statistiku, která reflektuje variabilitu ve směru jednotlivých vektorů. Na základě výše uvedeného tedy platí: Směrový rozptyl (circular variance) Sv se potom vypočte: Sv =l-OR/n kde n je počet vektorů. Svmůže nabývat hodnot 0 až 1. Je-li Sy=0, potom OR=n a všechny vektory mají stejný směr. Je-li Sy= 1, potom OR=0, všechny vektory mají opačný směr a výsledný vektor je bod. 6.5 Úvod do statistického popisu sítí Nebude probírána síťová analýza - ta vyžaduje speciální prostředí a nástroje (maticový počet) i speciálně upravená vstupní data. Základní pojmy používané v síťové analýze: nódy a hrany (spoje), jejich počet také charakterizuje síť. Ke křížení dvou a více hran dochází pouze ve vrcholu (planar graph topology) Deskriptory sítě lze rozdělit do dvou skupin: 1. Deskriptory sítě j ako celku 2. Deskriptory relací jednotlivých segmentů sítě. 6.6 Konektivita a matice konektivity Matice konektivity (tab. 3.1) shrnuje informaci o tom, které segmenty sítě spolu souvisí (jsou bezprostředně spojeny). Lze však charakterizovat i úroveň konektivity sítě jako celku. Pro fixní počet vrcholů má síť s větším počtem spojů lepší konektivitu. Dále existuje minimální počet spojů, který zajišťuje spojení všech vrcholů. Bude-li v - počet vrcholů sítě, e - počet hran sítě potom: e = v — 1 mm 79 Minimálně propojená síť (Minimally conneted network - MCN) - odstraníme-li jakoukoliv jednu hranu, síť se rozpadne na dva subsystémy. Podobně lze pro daný počet vrcholů vytvořit maximální počet hran, které spojují všechny vrcholy. Tedy maximální počet hran v síti o v vrcholech: =3(v-2) Jednoduchou charakteristikou konektivity sítě je Gamma index (y) - je definován jako poměr aktuálního a maximálního počtu vrcholů sítě. e y =- e max Další jednoduchou charakteristikou konektivity sítě je počet okruhů. Výskyt okruhů v síti značí možnost dostat se z jednoho místa do jiného alternativními cestami. Síť s minimální konektivitou nemá žádný okruh. Počet okruhů lze zjistit tak, že od aktuálního počtu hran v síti odečteme počet hran potřebný pro minimálně propojenou síť (MCN), tedy e-(v-l) nebo e-v+1. Obdobně pro daný počet vrcholů je maximální počet okruhů roven 2v-5. S oběma uvedenými počty okruhů lze vytvořit poměr aktuálního počtu k počtu maximálnímu - tedy tzv. alfa index: e-v + 1 a =- 2v-5 Pomocí alfa indexu můžeme snadno porovnat dvě sítě. 6.7 Dostupnost sítě (Acccessibility) Jedná se o charakteristiku jednotlivých vrcholů či hran sítě. Popisuje jejich dostupnost v rámci sítě. Další text se týká dostupnosti hran sítě, obdobné vztahy lze definovat i pro vrcholy. Jednoduchým ukazatelem dostupnosti hrany v rámci sítě je, s kolika jinými hranami daná linie přímo souvisí. Tuto informaci lze vyčíst z binární matice konektivity, pokud tuto doplníme např. řádkovým součtem. Tabulka 6.2 Matice konektivity a dostupnost hran v rámci sítě ID 1 2 3 4 5 6 7 8 9 10 SUMA 1 0 1 1 0 0 0 1 0 0 1 4 2 1 0 1 0 0 0 0 0 0 0 2 3 1 1 0 1 1 0 0 0 0 0 4 80 40010110000 3 50011010000 3 60001100000 2 71000000111 4 80000001010 2 90000001100 2 1ol 1000001000 2 Tabulka 6.3 Charakteristiky dostupnosti sítě (viz. obr 3.5, 3.6, 3.7) počet kroků k počet dosažení celkový počet přímých nejvzdálenějšího přímých a ID spojů místa nepřímých spojů 14 3 15 2 2 3 19 3 4 3 16 4 3 4 21 5 3 4 21 6 2 5 28 7 4 4 18 8 2 5 25 9 2 5 25 10 2 4 20 81 Obr. 6.5 Dostupnost jednotlivých segmentů sítě charakterizovaná počtem přímých spojů Uvedená charakteristika však může být zavádějící, protože nebere v úvahu relativní (topologickou) polohu hrany v rámci sítě. Hrana může mít i pouze jeden či dva spoje, přesto může být snadno dostupná, protože se nachází uprostřed sítě (a naopak). Relativní pozici každé hrany v rámci sítě lze zjistit např. pomocí počtu hran, kterými se lze z daného spoje dostat do nej vzdálenějšího místa sítě. počet kroků do nevzdálenějšího místa /V5 Obr. 6.6 Dostupnost jednotlivých segmentů sítě charakterizovaná počtem kroků nutných k dosažení nej vzdálenějšího místa sítě. Diametr (poloměr) sítě - je to jedna (1) plus největší počet hran nutných k dosažení nej vzdálenější ho místa v síti. Kvalitu spojení dvou hran (vrcholů) definuje počet hran mezi nimi. Spojení mohou být přímá a nepřímá. Tedy počet přímých a nepřímých spojů, které jsou třeba, aby byla daná hrana spojena se všemi hranami ostatními. Nepřímé spoje lze vážit počtem kroků. Zřejmě platí, že čím větší je celkový potřebný počet spojů, tím hůře dostupná j e daná hrana. Celkový počet spojů (přímých i nepřímých) je mírou dostupnosti. 82 Obr. 6.7 Dostupnost jednotlivých segmentů sítě charakterizovaná počtem přímých a nepřímých spojů nutných k dosažení jakéhokoliv místa v síti 7. Prostorové uspořádání ploch Využití prostorové statistiky k popisu měr úrovně a variability geografických jevů spojených s plochami (polygony) má v řadě geografických disciplín dlouhou tradici (demografie, krajinná ekologie apod.). Studium prostorových vztahů může být zaměřeno na následující typy úloh: 1. porovnání prostorového uspořádání studovaného jevu s uspořádáním teoretickým (shlukovým, pravidelným či náhodným) 2. typologie prostorového uspořádání jevů (bez územní souvislosti) 3. regionalizace - seskupování jednotek (polygonů) do vyšších územně souvisejících celků 4. interpolace a vyhlazování areálových dat 7.1 Míry prostorového uspořádání ploch Prostorová autokorelace- hodnoty atributů ploch spolu korelují v závislosti na jejich vzájemné poloze. To je v důsledku podobných přirozených (přírodních) podmínek (např. produkce zemědělských podniků) či v důsledku přirozené spojitosti jevů. U prostorově autokorelovaných dat nejsou hodnoty atributů v prostoru náhodné, ale prostorově závislé. Tato vazba (autokorelace) může být pozitivní (shlukové uspořádání -sousední objekty mají podobné hodnoty) či negativní (u pravidelného uspořádání). V případě náhodného uspořádání - slabá či žádná prostorová autokorelace. Také v případě prostorové autokorelace lze měřit její sílu. 83 Obr. 7.1 Příklad pozitivní prostorové autokorelace (shlukové uspořádání - vlevo) a negativní prostorové autokorelace (disperzní uspořádání - vpravo) Prostorová autokorelace je významným ukazatelem k hodnocení dynamiky a časových změn v prostorovém uspořádání objektů a pro predikce. Další význam prostorové autokorelace spočívá ve skutečnosti, že řada statistických ukazatelů (např. regresní modely) požaduje splnění předpokladu náhodnosti výběru objektů a jejich vzájemné nezávislosti. Míry prostorové autokorelace tak mohou potvrdit či vyvrátit splnění uvedených předpokladů. 7.2 Matice prostorových vah (Spatial weights matrices) Prostorová autokorelace měří stupeň podobnosti atributů mezi danou plochou a plochami sousedními. Nejprve proto musí být vztahy sousedství jistým způsobem kvantifikovány. Máme plochu s n prostorovými jednotkami. Potom můžeme definovat n x n párů sousedství -maticí typu n x n. Každá prostorová jednotka je prezentována jedním řádkem a sloupcem. Každá hodnota v matici prezentuje prostorový vztah mezi jednotkami prezentovanými daným řádkem a sloupcem v matici. Buňky matice mohu nabývat různých hodnot v závislosti na způsobu definování sousedství (např. binární matice s 0 a 1 podle toho, zda jednotky spolu přímo sousedí či nikoliv, nebo - buňky nesou vzdálenost mezi centroidy obou jednotek. Protože hodnoty v buňkách představují váhy při výpočtu prostorové autokorelace, potom se sestavené matice označují jako matice prostorových vah). 7.2.1 Způsoby definování sousedství Označují se podle pohybu šachových figur (Rook's case - věž, Queen's case - Dáma) - viz. Obr. 7.2 Bezprostřední sousedé (se společnou hranicí, i jedním bodem v případě Queens case) jsou sousedé prvního řádu. Analogicky lze definovat sousedy vyšších řádů. A % C D I ~ X I F i G H A C ■1 V/ _ X " / i x i F i G 4 H Obr. 7.2 Způsoby definování sousedství 84 Vedle sousedství je další běžně užívanou mírou prostorové relace objektů jejich vzdálenost. Intenzita vztahu dvou vzdálených jednotek bude obecně menší než intenzita vztahu jednotek blízkých. Tato vzdálenost může být arbitrárne určena (na základě zkušenosti či povahy studovaného problému: např. k danému domu jsou sousedé definováni jako domy do vzdálenosti 1 km, výsledek potom ze vyjádřit v binární podobě). 7.3 Binární matice konektivity (BCM - binary connectivity matrix) Analogicky jako v případě linií - binární, čtvercová symetrická matice C s prvky Cy, 1 -sousedí, 0 - ne) /d Bitie venko Bitia měsfo Hodcnvn Břeclav 0.0000 1.0000 1.0000 1.0000 0.0000 1.0000 1.0000 Blansko TľÓOOO' g.'üÖüÖ .............iTöoö'ö" i'.öüöü OľOGOO .............Ö'.'ÖÖÖÖ' clüööo' Vyškov TľÓOOO' .............i'.'öööö' .............ÖIÖÖÖ" ö.öööö iTög'üö .............Ö'.'ÖÖÖÖ' .............r'ö'ö'ö'ö" Brno-město iTÖOÖ'Ö" 'i'.'öööö' .............ö.öööö' o ciücici .............ö'.'öööö' Hodonín öTö'öoö" ö'.'öööö' i"öoöö" .............ö'.'öööö' □ ciücici .............ö'.'öööö' .............t'ö'ö'öö' Znojmo iTö'ö'ö'ö" 0.0000 clöciöö .............ö'.'öööö' oľcicicici .............ö'.'öööö' .............t'ö'ö'öö' Břeclav iTö'ö'ö'ö" 0.0000 1 o.'öooo' ''iTöoö'ö' i'.'öööö' o.o'ooo 7.3.1 Binární matice sousedství Vlastnosti BCM: • Prvky na hlavní diagonále mají hodnoty 0 • Matice j e symetrická - redundance uložené informace • Suma v řádku nese informaci o počtu sousedů dané jednotky • Pro větší počet prostorových j ednotek obsahuj e velké množství nul a j e tedy paměťově náročná Vhodnější způsob zaznamenání vztahů sousedství je uchovávání ID či názvu sousedů pro každou plochu, tedy např.: Polygon Sousedi Soused2 ... Brno-město Brno-venkov Blansko Blansko Brno-venkov Vyškov Brno-město 7.3.2 Stochastická matice či matice se standardizovanými řádkovými vahami (RSWM) Zaznamenání sousedství v binární podobě není v řadě případů výhodné - váhy jsou stejné bez ohledu na počet sousedů. Vhodnějším způsobem je nahrazení jedniček vahou wij , vypočtenou jako poměr mezi hodnotu cij a sumou v řádku - tj. počtem sousedů. Tedy má-li jednotka 4 sousedy, bude její váha rovna 0,25 - tak dostaneme z matice C matici W, 85 označovanou jako matici se standardizovanými řádkovými vahami. Stejně jako matice C má i W na hlavní diagonále nuly, není vak již symetrická. Bma_ ¥&rski\ Blanska Vysl:av Brna měsfa Hadotxn Ztiajnia Brachy Brno^yenkgy 0.0000 0.2000: 0.2000i 0.2000 0.0000 0.2000 0.2000 Blansko 0.3333! 0.0000! 0.3333 i 0.3333 0.0000 0.0000 0.0000 Vyškov 0.2500! 0.2500! 0.0000j 0.0000 \ 0.2500 0.0000 0.2500 Brno-město 0.5000 í 0.5000! 0.0000 i 0.0000 [ 0.0000 0.0000 0.0000 Hodonín 0.0000 i 0.0000! 0.5000! 0.0000 [ 0.0000 0.0000 0.5000 Znojmo 0.5000 i 0.0000! 0.0000! 0.0000 [ 0.0000 0.0000 0.5000 Břeclav 0.2500! 0.0000! 0.2500! 0.0000 0.2500 0.2500 0.0000 Obr. 7.3 Matice se standardizovanými řádkovými vahami 7.4 Vzdálenosti centroidů Vztahy prostorové závislosti lze charakterizovat také vzdáleností jednotek (viz. první zákon geografie - Tobler, 1970: Všechny objekty spolu souvisí, ale blízké objekty spolu souvisejí více). Tedy vzdálenost je vhodnou váhou pro definování prostorových vztahů. Existuje několik způsobů definování vzdálenosti dvou polygonů, např. vzdálenost centroidů. Existuje několik způsobů určení centroidů pro daný polygon. V závislosti na tvaru polygonu nemusí jeho centroid ležet uvnitř něho. Jsou-li jako váhy použity vzdálenosti (zde vzdálenosti centroidů), matice se označuje D s prvky dij . Váhy jsou potom definovány jako převrácená hodnota vzdálenosti: 1 7 d V řadě případů síla vztahu mezi dvěma jednotkami klesá rychleji než se zvětšuje jejich vzdálenost, proto se váhy definují jako. 1 147.. = - 7 d2 7.4.1 Nejbližší vzdálenosti Na místo vzdáleností centroidů jsou použity vzdálenosti dvou nejbližších částí dvou polygonů. Takto definované váhy jsou výhodné pro charakterizování prostorových kontaktů či difúze. U takto sestavené matice buňky s nulami mimo hlavní diagonálu (sousedé) odpovídají buňkám s jedničkami v binární matici sousedství. 86 Brna město Hadotxn i Brno-venkgy 0.0000 0.0000 0.0000 0.0000 6.3679 0.0000 0.0000 ' Blansko 0.0000 0.0000 0.0000 0.0000 23.0282 29.5297 24.4276; Vyškov \ 0.0000 0.0000 0.0000 17893 0.0000 23.7376 0.0000; Brno-město ..............aoooo' 0.0000 17893 0.0000 157463' 14*2933 161121 Hodonín 13679* 210282 aoooo 157463 CL0000 '30.5051 a'oooo] Znojmo aoooo" '215297 217376' 14*2933 315051' 0.0000 a'oooo] Břeclav *0.0000 24.4276 110000 16112 '"' 0.0000 '"' 0.0000 '"' 0.0000 Obr. 7.4 Matice vzdáleností mezi nejbližšími částmi polygonů 7.5 Míry prostorové autokorelace Výše uvedené matice slouží k definování měr prostorové autokorelace (SA). Míry SA mohou být vztaženy k poli bodů (viz. výše) či ploch. V případě ploch lze zpracovávat data nominální (JCS - joint count statistics - Statistika charakteru sousedství), intervalová i poměrová (Moranův index I, Gearyho poměr C, G-statistika) Uvedené míry lze označit jako globální míry prostorové autokorelace (asociace). Tedy jedna hodnota je vypočtena pro celou studovanou oblast. Avšak také prostorová autokorelace se může měnit v rámci studované oblasti - k deskripci prostorové heterogenity prostorové autokorelace lze využít lokálních měr - Local Indicator of Saptial Association (LISA) a lokálmí verze G-statistiky (local G-statistics). Ke grafickým prostředkům hodnotícím prostorovou autolorelaci patří Moranův scatterplot diagram. Základní notace používaná v následujícím popisu indexů prostorové autokorelace Wjj - obecně buňka matice vah Wpro řádek i a sloupec j. (nejen matice stochastické - viz. výše) Sumace vah daného řádku i přes všechny sloupce (řádková suma): i Sumace vah daného sloupce j přes všechny řádky (sloupcová suma): i Sumace všech buněk matice vah: i j Pro testování významnosti indexů prostorové autokorelace lze váhy v jednotlivých výrazech sumarizovat do následujících výrazů: 87 suH=^Z2>*i+wii)2 1 i j f Y SUMi - suma přes váhy. Jsou-li váhy binární a matice symetrická, potom (v^ + wj;)2 = 4 SUMi je tedy čtyřnásobek celkového počtu spojů (společných hranie) v celé studované ploše. Hodnota SUM2 je založena na sumování vah každé plošné jednotky v obou směrech (wij i wjí). Výsledná hodnota je potom získána jejich součtem, umocněním a sumací pro všechny jednotky studované oblasti. Nechť n je počet plošných jednotek ve studované oblasti. Existují-li dvě skupiny jednotek definovaných atributy s hodnotami x a y, potom výrazy nx a ny značí počet jednotek v jednotlivých skupinách. Podobně: n(x) =n*(n-1) * (n - 2) * (n - 3) *...*(«- x +1) kde n > x Například, bude-li n=5, potom n(3) = n(n-l)(n-2) = 5x4x3 a n(1) = n Jestliže xi je hodnota atributu pro plochu i, můžeme definovat nový parametr m,, založený na hodnotách xi: i=l kde j = 1,2,3,4. Potom, jestliže j= 1, m, je suma xi pro všechna i. Jestliže j=2, m, bude suma všech čtverců xi. 7.5.1 Statistika charakteru sousedství - Joint count statistics (JCS) Touto metodou lze zjistit, zda uspořádání ploch, které mohou nabývat binárních hodnot vykazuje prvky náhodnosti. Tedy zda existuje pozitivní (clustered pattern) či negativní (random pattern) prostorová autokorelace. 88 Obr. 7.5 Statistika četnosti spojů (JCS) Podstata metody - jednoduchý příklad: Máme mapu se dvěma kategoriemi landuse: U - zástavba, R - volná krajina. Potom mohou existovat čtyři typy sousedských vztahů: UU, RR, UR, RU. V případě čistě náhodného uspořádání se bude každá kombinace vyskytovat v 25% případů. Dvojice ploch s odlišným atributem se budou vyskytovat v 50 % případů. Pokud UR + RU < 50%, potom výskyt dvojic ploch se stejným atributem UU a RR bude vyšší než 50% - což je případ pozitivní prostorové autokorelace. V případě 50 na 50 - uspořádání je náhodné a pokud UR + RU > 50%, pak se jedná o negativní S A, kdy dominují hranice nepodobných ploch. Mapu (obr. 1) s pěti plochami můžeme prezentovat také grafem s vrcholy a spoji, zaznamenávajícími druh povrchu a také bezprostřední sousedství jednotlivých ploch s plochami jinými, jak je patrné z obr. 4.4 Obr. 7.6 Grafická prezentace druhů spojů Sestavíme matici sousedství pro jednotlivé plochy. V této matici nula značí, že obě plochy spolu bezprostředně nesousedí, 1 naopak. Zároveň je barvou buňky v matici naznačeno, o jaký typ spoje se jedná (Obr. 7.7). □□□□□ □□□□□□ □□□□□□ □□□□□□ Obr. 7.7 Binární matice sousedství pro nominální data 89 Pořadí řádků a sloupců v uvedené matici je určeno abecedním pořadím identifikátorů ploch. Nic nebrání sestavit matici v jiném pořadí řádků a sloupců - například podle typu povrchu -(viz Obr. 7.8). □□□□□□ □□□□□□ mmra Obr. 7.8 Binární matice sousedství uspořádaná podle hodnot atributů Obě matice jsou symetrické, ve druhém případě navíc je možné jednoduše popsat prostorovou autokorelaci pomocí čtyř sub-matic. Z matice lze zjistit, že 14 buněk obsahuje jedničku, která značí výskyt hrany (14 párů sousedství). Dále platí, že jednotlivé typy sousedství se na mapě vyskytují s těmito četnostmi: • UU=2 • UR=5 • RU=5 • RR=2 Z toho plyne, že RU + UR > 14/2 , tedy naše mapa vykazuje negativní autokorelaci, nepodobné plochy (s odlišným typem povrchu) se shlukují. Uvedený koncept lze dále rozšířit využitím počtu pravděpodobnosti a statistických testů. Ty nám umožní testovat statistickou významnost prostorového uspořádání ploch v mapě. V dalším výkladu jsou používány dvě hodnoty atributů B - black, černá, W - white, bílá. Tedy bude-li prostorové uspořádání indikovat uspořádání do shluků, potom můžeme předpokládat více hranic typu BB či WW než BW nebo WB - tedy pozitivní prostorovou autokorelaci. JCS tedy nejprve určuje počet jednotlivých druhů spojů s cílem testovat četnost jejich výskytu. Pro plochu s malým počtem polygonů lze počty jednotlivých spojů zjistit manuálně, pro velký počet ploch je nutné využití metod matematické statistiky. Obecné kroky výpočtu jsou následující: Nechť xi= 1 jestliže polygon i je černý a xi=0 jestliže polygon i je bílý. Potom pro BB spoje bude: 0BB = ^^^(w^Xj) Pro WW spoje bude platit: Cw = £. [wi} (1 - jq)(1 - Xj)] Pro BW nebo WB spoje bude platit: 0BW =~zZ zZ K ~ xj )2] Uvedené vzorce představují výrazy pro pozorované (O - observed) počty spojů popisující dané uspořádání. Vysoké hodnoty Obb či Oww či obou indikují pozitivní prostorovou autokorelaci (slukování). Pozorované počty spojů však musíme porovnat s náhodným uspořádáním a musíme testovat, zda eventuelní zvýšené počty Obb či Oww nejsou výsledkem pouhé náhody, zda jsou či nejsou statisticky významné. Budeme tedy pracovat s počtem pravděpodobnosti. Způsob určení pravděpodobnosti výskytu B a W polygonů však může významně ovlivnit výsledek analýzy. Hodnoty atributů mohou byt jednotlivým polygonům přiřazeny na základě předpokladu normality či náhodnosti (viz. prostorová analýza bodů) Předpoklad normality: (NORMALITY - FREE - SAMPLING) - pravděpodobnost, že se jedná o polygon B či W je založena na teorii či na trendu hodnot atributů odvozeném z větší oblasti. Pravděpodobnost, že polygon má B či W není ovlivněna celkovým počtem B či W polygonů v oblasti. Předpoklad náhodnosti: (RANDOMIZATION - NONFREE - SAMPLING) -pravděpodobnost, že polygon bude mít B či W je omezena či závisí na celkovém počtu B či W polygonů. Příklad: Plocha obsahující sedm polygonů: Předpoklad náhodnosti - může existovat různá konfigurace 4 „černých" a 3 „bílých" ploch. Předpoklad normality - může existovat různá konfigurace jakéhokoliv (0 až 7) počtu „černých" a „bílých" ploch. U metody JCS bychom neměli pracovat s předpokladem normality v případě, že informace získané z teorie, zkušenosti či z trendové funkce z širšího okolí jsou nespolehlivé. Náhodné vzorkování totiž vyžaduje méně rigorózní podmínky použití. 7.6 Normální vzorkování V obou výše komentovaných případech je nutné vedle pozorovaných (O) počtů jednotlivých typů spojů či hranic (joint) zjistit počty očekávané (E) a také jejich směrodatné odchylky. Očekávané počty odrážejí efekt náhodnosti či nevýznamné prostorové autokorelace jakéhokoliv typu (pozitivní či negativní). Tedy zjistí se diference mezi pozorovanými a očekávanými četnostmi spojů. Tyto diference j sou následně standardizovány hodnotami příslušných směrodatných odchylek a získáme tak standardizovaná skóre. Z hodnot těchto skóre můžeme rozhodnout, zdaje ve studované oblasti významná pozitivní či negativní prostorová autokorelace v uspořádání polygonů podle hodnot atributu. Jinými slovy, je nutné provést tři typy porovnání. Dále je prezentován případ pouze pro testování negativní prostorové autokorelace. Pro případ normálního vzorkování j sou vztahy pro očekávané četnosti jednotlivých druhů spojů (joint) (EBb, Eww, EBw) následující: 91 EBB = |WP2 Eww = |wq2 EwB=Wpq p - pravděpodobnost, že plocha bude B (černá) q - pravděpodobnost, že plocha bude W (bílá) Pravděpodobnosti p, q musí dávat 100% nebo (p + q = 1). Pokud není k dispozici jiná informace, potom p = nB/n, jsou však i jiné způsoby určení p. Pokud je použitá prostorová matice vah binární, lze výrazy pro očekávané počty typů spojů zjednodušit: kde J značí celkový počet spojů ve studované oblasti. K testování statistické významnosti zjištěného prostorového uspořádání lze využít Z-testu. K němu je zapotřebí zjistit směrodatné odchylky očekávaných počtů spojů. Směrodatné odchylky se vypočtou v závislosti na použité váhové matici následovně: Pro stochastickou matici vah: ebb = JP2 Eww = Jq2 EBW = 2Jpq Pro binární matici vah: bb Vp2J + p3K-p4(j +K) ww Vq'J +q3K-q4(j +K) kde a je směrodatná odchylka počtu příslušných spojů Si, S2, J, p, q byly definovány výše 92 Hodnota n v tomto výrazu značí celkový počet polygonů a Li je počet spojů mezi polygonem i a jeho sousedy. Obecný postup testování (na príkladu negativní prostorové autokorelace (BW spoje) při použití binárni matice): Pro výpočet očekávaných potřebujeme znát hodnoty pravděpodobností p, q. Rozhodneme se pro určité pravidlo definující sousedství (rook, queen). Dále určíme J (počet spojů) - zjistíme sumováním všech členů binární matice vah a dělíme dvěma. Odhad správných hodnot p a q -ze zkušenosti, z teorie (např. mortalita v určitém regionu - použijeme údaje o mortalitě celého státu. Potom určíme hodnotu výrazu L(L-1) pro každý polygon a provedeme sumaci pro celou oblast. Potom určíme hodnoty EBwa obw- Máme-li k dispozici pozorované počty spojů (0Bw), potom můžeme vyjádřit hodnotu z-skóre: O - E 2 — wbw ^bw ^bw Podle pravděpodobnosti rozdělení hodnot Z-skóre platí, že jakákoliv hodnota Z ležící mimo interval (-1,96; -1,96) má pravděpodobnost výskytu menší něž 5 případů ze 100 (a=0,05). (a) Clustered: 4 BWjoins (c) Dispersed: 8 BWpins (b) Random: 6 BW joins (d) Number of joins Obr. 7.9 Příklady prostorového uspořádání černých a bílých polygonů v rámci studované oblasti (a, b, c) a počty sousedů jednotlivých ploch (d) PŘÍKLAD: Na obrázku (Obr. 7.9) je oblast obsahující 7 polygonů. Naším cílem je metodou JCS určit, zda v této oblasti existuje statisticky významná negativní prostorová autokorelace ve výskytu „černých" (B) a 93 „bílých" (W) ploch. Jako vah využijeme prvků binární matice. Podle výše uvedených vzorců musíme vyčíslit hodnoty 0Bw, EBw obw, 1) Spočteme celkový počet všech spojů ve studované oblasti, tedy hodnota J= 11. 2) Určíme způsob definice sousedství - v tomto případě za sousedy považujeme pouze polygony, které spolu sousedí hranou (rook's case). 3) Určíme hodnoty pravděpodobností p, q výskytu „černé" či „bílé" plochy. V tomto případě předpokládáme, že p=0,3 a q=0,7. 4) Z obr. d určíme pomocí následující tabulky hodnotu ~S\l(L -1) Oblast L L-l L(L-l) A 3 2 6 B 2 1 2 C 3 2 6 D 5 4 20 E 3 2 6 F 3 2 6 G 3 2 6 E 22 52 5) Vyčíslíme hodnoty , EBw 0bw: EBW =2Jpq = 2*11*0,3*0,7 =4,62 <7bw —2,1 6) Pro jednotlivé varianty na obrázku a, b, c jsou hodnoty pozorovaných počtů spojů (Obw) Obw =4,6 resp 8 7) Pro konfigurace „černých" a „bílých" poch uvedené na obrázku vyjádříme hodnotu z-skóre: a) Z=±^= = -0,29 2,1 b) Z=^ = 0,65 C)Z=Í^= = 1.61 2,1 8) Interpretace: Žádná z hodnot Z-skóre nepřesahuje prahovou hodnotu ±1,96 a tedy uvedená uspořádání nevykazují statisticky významnou negativní prostorovou autokorelaci na hladině významnosti a=0,05. 94 7.7 Náhodné vzorkování V tomto případě závisí pravděpodobnost, zdaje polygon bílý nebo černý, na celkovém počtu černých polygonů a počtu bílých polygonů ve studovaném území. Obrázek 4.7. uvádí tři typy prostorového uspořádání sedmi polygonů ve studované oblasti. Protože ve všech třech případech j sou počty B a W polygonů stejné (jsou jen jinak uspořádané) hodnoty pravděpodobnosti budou: p=3/7 a q=4/7. Dále se vypočtou hodnoty očekávaných počtů spojů a jejich směrodatné odchylky. Výpočetní vzorce jsou jiné než v případě normálního vzorkování (viz. Lee, Wong, 2000, str. 154 - 155). Postup výpočtu je však analogický výše uvedenému příkladu. V.ttlttical AiutlytH mth Ar< Vww E*> E« y*w Itww SpiMAutoamMon Qiaehct ttndow tM> |01MMb|(žtmál/l,iMTJ-.l Se* 1:12.898.834 j. Otxn Oarti ll> Sort AAJorl!>4 88 Jcrts-2 «8 JwHi ■ 7 EnpecWjAAJcrti-4 2449 Enpected 88 Jortt - 2 38778 E»tecletä \ Sctmnvt Portage Asltfahit/a Lake \ 1 Geauga 0.0000 25.1508 26.7057 32.7509 25.0389 26.5899 12.6265. Cuyahoga 25 i508 ............0.0000 47 8151 3Í 6155 50.80G4 28 2214' Trumbull 26 7057 47.8151" ' 0 0000 ' 41.8561 24 4759 '29.5633 36.7535! Summit 32 7509 '23.4834 41 8561 b.boob 17 8031' '58.0869' 42.7375: Portage 25 0389 '31.6155 24 4759 ' 17.8031 0 0000 '45.5341 37 4962 ■ Ashtabula 26 5899 '50.8064' 29 5633 ' 58.0869 45 5341' ............o.'obbo 24.7490: Lake .2.. 28.2214 ' 42.7375 37.4962 b.o'ooo: «l T Obr. 7.12 Výchozí matice vzdáleností centroidů w. distmatriH.dbf Id I Geaitga \ Cciyaliaga \ Tfwihiá \ Swwxt \ Poftage \ Asítfahida \ Lake Geauga Cuyahoga Trumbull Summit Portage Ashtabula Lake 0.0000 J TobooT 1 00001 0 0000 j 1 oooo! 1 oooo! ' i oooo! 1.0000! ' o!oooo ľ '""oiooooi TobooT """"biooooT """ oioooo'ľ TobooT 1.0000] TobooT "abbbb'T TobooT "TobooT Tbbbb'T "abbbb'T o.oooo: TboboT □"□OOO; oToooä ľ Trabl □ToioiooT □Toobo ľ 1.0000] —I' "iľobbbT "TobooT 'aoooo'T "abbbbj '"abbbb! 1.0000: TobooT TooooT "b^bbbbT "oioboo'T '"oioboo'T '"í"bbbbT 1.0000 j ■■■——■j □.□□bo] '"obbbbl 'abäoo'1 "iľbbbbl '"obbbbl Obr. 7.13 Matice sousedství vypočtená pro d=30 z matice na obr. 5.2 MM G-Statistics = 0.555756 The Expected G = 0.52381 The Variance of G = 0.00856308 Z-Value of G = 0.345226 ■i Obr. 7.14 Výsledky výpočtu obecné G- statistiky pro vstupní data na obrázku 5.1 pň použití matice vzdáleností centroidů a hodnotě definující vzdálenost d=30 mil. Vypočtená hodnota G(d) vykazuje mírnou úroveň prostorové asociace, podle hodnoty z-skóre však výsledek není statisticky významný. Jinými slovy - dané uspořádání průměrného přijmu v sedmi státech Ohia je spíše výsledkem náhody než určitého systematického procesu. 7.9.1 Lokální statistiky prostorové autokorelace Všechny tři uvedené indexy jsou příkladem indexů globálních. Jsou sumární hodnotou prostorové autokorelace pro celou zpracovávanou oblast. Je však pravděpodobné, že hodnoty prostorové autokorelace se budou v různých sub-oblastech měnit. Navíc můžeme očekávat, že pozitivní autokorelaci lze nalézt v jednom sub-regionu a negativní v jiném. Proměnlivost prostorové autokorelace v rámci studované oblasti lze vyšetřovat výše uvedenými indexy modifikovanými pro detekování prostorové autokorelace v lokálním měřítku. LISA (Local Indicators of Spatial Association) Jedná se o lokální verze Moranova a Gearyho indexu. Ke zjištění úrovně prostorové autokorelace na lokální úrovni je nutné vypočítat hodnotu indexu pro každou plochu zpracovávaného území. Lokální Moranův index pro jednotku i je definován takto: t^z1 kde ži a Zj jsou odchylky od průměru nebo 100 7_U-x) (7 kde o j e směrodatná odchylka xi. Podobně jako v případě globálního Moranova indexu znamenají vysoké hodnoty kumulaci podobných hodnot atributů (vysokých či nízkých) v sousedních plochách, nízké hodnoty potom kumulaci odlišných hodnot atributů. Obecně hodnoty Wy mohou představovat po řadách standardizovanou matici vah, lze použít i jiných matic vah. Zjištěné hodnoty lokálního Moranova indexuje nutné porovnat s očekávanými hodnotami a testovat statistickou významnost jejich rozdílu pomocí z-skóre. Očekávané hodnoty při hypotéze náhodnosti: E[li] = -wi./(n-l) a hodnota rozptylu: Var[l\- (II- -m4 n-\ Iml) + 2w. i(kh) (2m4/m22 -n) wf (n-\)(n-2) (n-lf kde W2 W; v j j a výraz 2wíw=ZZwikwih k^i h^i Každá plocha ve zpracovávaném území má svoji I hodnotu a té přísluší hodnota očekávaná a také jistá hodnota rozptylu. Hodnoty I mohou být vynášeny do mapy v podobě kartogramu. Lokální verze Gearyho poměřuje definována následovně: s =Ewíj(zí-zj)2 Hodnoty rozdělení lokálního Gearyho indexu nemají tak vhodné vlastnosti jako v případě indexu Moranova. Jejich interpretace je však obdobná jako v případě globální verze indexu. Shlukování podobných hodnot atributů vede k nízkým hodnotám tohoto indexu a naopak. 101 Lokální G-statistika Měří asociaci hodnot atributů v ploše i a v plochách okolních definovaných vzdáleností d: G;(d)=——-pro i*j Obdobně jako v předchozích případech je nutné interpretovat hodnotu indexu pomocí, očekávaných hodnot, hodnot rozptylu a standardizovaných skóre. Očekávané hodnoty se vypočtou následovně: E(Gi)=\Vf/(n-l) kde Definice rozptylu: Var(Gi) = E(Gi2)-[E(Gi)]2 E(G;2) (Z^)2 "^(n-l-^)£.x? (n-l)(n-2) (n-l)(n-2) pro i* J Vysoká hodnota z-skóre je spojena s výskytem shluků podobných a vysokých hodnot indexu. Jestliže je shluk tvořen nízkými hodnotami, z-skóre bude nabývat velkých záporných hodnot. Hodnoty z-skóre kolem nuly indikují neexistenci zřejmého prostorového uspořádání hodnot atributů v plochách studovaného území. Příklad 3: Pro data z příkladu 1 byly vypočteny hodnoty lokálního Moranova indexu I (pro každý stát). Jako matice vah byla použita matice stochastická (Obr. 7.15). Výsledky jsou prezentovány ve formě kartogramu na následujících obrázcích (Obr. 7.16 a Obr. 7.17). 1 0 distmatrix.dbf -M *l 1 M \ Seauqa Qtyafjoija \ Tritnibtä Stjnwiit Portage \ Ax/tŕeĎďa Lake i Geauqa ! o 0000 0.1667] 0.1667 0.1667 0.1667 0.1667 0.1667 Cuyahoga C 2500 0 0000 ľ ooooo' 0.2500 ' 0 2500 'ooooo "o! 2500 Trumbull C .... ooooo ľ ooooo' 0.0000 ' 0 3333 '0 3333 " o! 0000 Summit C .... 0 3333i' OOOOO' 0.0000 ' 0 3333 'ooooo " o! 0000 Portage C 2500 0 2500 i' " 0 2500' 0.2500 ' 0 0000 'ooooo " o! 0000 Ashtabula 0 .... 0 0000 i' " 0 3333' 0.0000 ' 0 0000 'ooooo 03333 Lake j 0 .... 0 3333i' " OOOOO' 0.0000 ' 0 0000 '0 3333 " o! 0000 «1 Obr. 7.15 Stochastická matice vah k definování sousedství pro výpočet lokálního Moranova indexu I 102 Moran I -1.012 -1 012- -0 634 -0.634--0172 ■0.172--0 02 •0.02-0.051 Obr. 7.16 Kartogram hodnot lokálního Moranová indexu I Z-skóre ~~] -1.892 ■ 1.892 1.693 B 1.693 -0.267 | -0.267 - 0.176 ■ 0.176-0 484 Obr. 7.17 Kartogram hodnot z-skóre pro lokální Moranův index I Interpretace: Vysoké hodnoty indexu I mají ty státy, jejichž sousedé mají velmi podobné hodnoty studované charakteristiky. Podle z-skóre žádná z hodnot není statisticky významná a dané uspořádání průměrných příjmů v sedmi státech lze interpretovat jako náhodný proces. Obdobným způsobem lze vizualizovat a hodnotit výsledky analýzy založené na lokálním indexu C a lokální G-statistice. Moranovo korelační pole (Moran Scatterplot) Lokální statistiky vystihují prostorovou heterogenitu v jednotlivých částech studovaného území. Pomocí nich je tedy možné jistým způsobem identifikovat oblasti s neobvyklými hodnotami měr prostorové autokorelace, které lze označit jako oblasti s odlehlými hodnotami (outliers). Efektivním nástrojem pro takovouto diagnostiku území je Moranovo korelační pole založené na regresním počtu. Předpokládejme, že x značí vektor hodnot xi s odchylkami od průměru (x; - x) a dále W značí po řádcích standardizovanou matici vah. Potom můžeme sestavit regresní závislost hodnot Wx na x. Směrnice této regresní závislosti indikuje vzájemný vztah sousedních hodnot atributů. Tedy x = a + IWx 103 kde a značí vektor koeficientů - (intercept). Hodnota I je regresní koeficient reprezentující směrnici a také hodnotou Moranova globálního indexu I. Vynesení regresní závislosti Wx na x umožňuje identifikovat odlehlé hodnoty. Pokud budou mít všechna pozorování podobné hodnoty prostorové autokorelace, v korelačním poli budou body blízko regresní čáry. Naopak pokud některá pozorování budou ukazovat lokálně výrazně vysoké či nízké hodnoty prostorové autokorelace ve vztahu k jejich sousedům, tato pozorování budou v grafu tvořit body výrazně nad či pod regresní čarou. Regresní čára vyjadřuje obecný trend hodnot prostorové autokorelace v celém zpracovávaném území a parametr její směrnice je index I. Příklad 4: Hodnota Moranova indexu (viz. Příklad 1) indikuje slabou negativní prostorovou autokorelaci (státy s vysokou hodnotou studovaného atributu j sou blízko států s nízkými hodnotami). Moran ScatterPIot for Medhinc89 : R-square = 0 816821 0.8 0.6 0.4 0.2 0 -0.2 -0.4 □ a = 0.261394, b =-0.305B48 .5 -1 -0.5 0 0.5 1 1.5 2 Obr. 7.18 Výsledek regresní analýzy a Moranovo korelační pole (Moran Scatterplot) pro průměrný příjem sedmi států Ohia (příklad 1). Parametr b představuje hodnotu Moranova indexu I Z grafu je patrné že příjem (x) je nepřímo úměrný vážené hodnotě příjmu (Wx). Množinou bodů lze proložit přímku. Body, které se výrazně odchylují od přímky představují „outliers" -představují oblasti s výrazně odlišnými hodnotami prostorové autokorelace. Interpelace s ohledem na polohu bodů v jednotlivých kvadrantech • high-high,low-low (2. nebo 3. kvadrant) = spatial clusters • high-low,low-high (1. nebo 4. kvadrant) = spatial outliers 104