SOUHRN

Od prosince 2012 je pro potřeby operativní hydrologie v Českém hydrometeorologickém ústavu (ČHMÚ) pravidelně pro každou zimní sezonu určována poloha nulové izochiony (sněhové čáry). Důvodem je odhad zásob vody ve sněhové pokrývce, s nimiž čeští hydrologové musejí nutně pracovat, pokud chtějí, aby jejich předpovědní modely poskytovaly relevantní výsledky. Pro lepší představu o aktuální prostorové distribuci sněhové pokrývky v Česku je informace o nulové izochioně odvozována ze snímků MODIS pořízených družicí Terra. Získaná časová řada reprezentuje již poměrně dlouhé období (nyní až do května 2021), a tak se nabídla možnost analýzy prostorové a časové dynamiky nulové izochiony v Česku. V této studii byla informace o izochioně rozdělena do 27 geomorfologických oblastí, přičemž zimní sezona byla začleněna také do období akumulace a tání sněhu. Snahou bylo zjistit, jaké jsou rozdíly mezi oblastmi a jednotlivými obdobími a jaké jsou vazby dynamiky nulové izochiony na vybrané faktory odvozené z dalších geografických dat, jako je digitální model reliéfu apod. Data o izochioně byla z nejrůznějších důvodů neúplná a nevyhovovala nasazení modelů, které vyžadují pravidelné rozestupy v čase. Proto bylo přistoupeno k odhadu chybějících denních hodnot tak, aby byly pravidelně pokryty zimní sezony vždy od listopadu do května. To bylo provedeno vhodnou modifikací EM algoritmu, jež zohledňuje jak strukturu časové řady, tak prostorové vazby. Následně byla aplikována korelační a regresní analýza, při níž bylo hlavním cílem zjistit, do jaké míry má vliv příslušnost ke geomorfologické oblasti (s jejími vybranými atributy) a zda dochází k signifikantním meziročním změnám.

ÚVOD

Nulová izochiona je již více než deset let praktickým pomocníkem při výpočtech zásob vody ve sněhu v rámci pravidelných hydroprognózních analýz pod hlavičkou ČHMÚ. Hlavním přínosem určování této izochiony je vymezení prostoru, kde lze předpokládat výskyt sněhové pokrývky a kde naopak sněhová pokrývka bude chybět, respektive kde lze počítat s vodní hodnotou sněhu a kde nikoli. Definování takové hranice v rámci homogenních regionů Česka umožňuje lépe interpolovat hodnoty výšky sněhové pokrývky zaznamenané ve staniční síti ČHMÚ. Bližší popis hydroprognózních analýz, včetně aplikace nulové izochiony při výpočtech, nabízí [1, 2]. Obsah tohoto článku navazuje na práci [3], kde byly shrnuty základní postupy extrakce nulové izochiony pomocí družicových snímků a informace o jejím rozložení v rámci geomorfologických oblastí Česka. Právě hledání korelací a dalších vazeb, které se dají ze sledované řady vypozorovat, jsou nosným tématem následujícího textu. Jednou ze stěžejních oblastí studia je kvantifikace míry závislosti variability průměrné nadmořské výšky nulové izochiony na terénních charakteristikách geomorfologických oblastí. Analyzována je jak dlouhodobá nadmořská výška nulové izochiony za celé zimní období vymezené měsíci listopadem až květnem (pro devět let 2013–2021), tak i v částech zimní sezony typických pro akumulaci sněhu, pro jeho tání a ve zbylé části této sezony (pro čtyři roky 2018–2021). Snahou bylo vysledovat míru působení faktorů souvisejících s konfigurací terénu v jednotlivých obdobích zimní sezony. Jelikož délka nasbírané časové řady je již poměrně dostačující (ve smyslu nečleněné zimní sezony), neméně důležitým úkolem bylo zjistit, jak se mění nadmořská výška nulové izochiony s časem, tj. v jednotlivých letech, a zda lze vypozorovat významný trend u některých geomorfologických oblastí. Při extrakci nadmořské výšky nulové izochiony a respektování její definice dle [4] bylo z důvodu nutné konzistence postupováno obdobně jako v [3]. Řešení nových úkolů bylo provedeno s využitím nejrůznějších statistických technik, přičemž kromě deskriptivní statistiky a metod pro doplňování chybějících hodnot hrála ústřední roli regresní analýza a výběr významných vysvětlujících proměnných.

DATA A METODIKA

Družicová data

Podstatná část dat analyzovaných v rámci tohoto projektu pochází z portálu Národního střediska sněhu a ledu (National Snow and Ice Data Center, dále NSIDC), které podporuje výzkum kryosféry, tedy sněhu, ledu, ledovců a zmrzlé půdy, ale i klimatických interakcí, jež v kryosféře probíhají. NSIDC spravuje a distribuuje vědecká data, vytváří nástroje pro přístup k datům, podporuje uživatele dat, provádí vědecký výzkum a vzdělává veřejnost o kryosféře. Jako platforma dat pocházejících z Národního úřadu pro letectví a vesmír (National Aeronautics and Space Administration, NASA) je zároveň certifikována coby tzv. CoreTrustSeal s osvědčením Pravidelný člen Světového datového systému, mezioborového orgánu Mezinárodní vědecké rady (International Science Council, ISC; dříve ICSU). Portál distribuuje data bezplatně celé vědecké komunitě již od roku 1976 v rozmanitých formátech známých v oblasti DPZ a ve velikostech od malých textových souborů po terabyty dat. Bližší informace o produktech, nástrojích i publikovaných výstupech lze dohledat na [5].

Konkrétní datová sada využitá pro účely analýzy nulové izochiony je označena jako MODIS/Terra Snow Cover 5-Min L2 Swath 500m, Version 61 a je včetně metadat dostupná z [6]. Název sady obsahuje základní popisné údaje o snímaných datech. Snímky jsou sbírány pomocí senzoru MODIS (Moderate Resolution Imaging Spectroradiometer), který je instalován na družici Terra. Terra je nadnárodní vědecký výzkumný satelit NASA na sluneční synchronní oběžné dráze kolem Země, jenž provádí simultánní měření zemské atmosféry, půdy a vody, aby přispěl k pochopení, jak se Země mění, a identifikoval důsledky pro život na Zemi [7]. Umístění senzoru MODIS na družici Terra ukazuje obr. 1.

Obr. 1. Družice Terra, vypuštěná 18. prosince 1999 (výška oběžné dráhy: 713 km; rychlost na oběžné dráze: 7 503 km·s-1; maximální rychlost: 27 010 km·h-1) a poloha senzoru MODIS (zdroj: [8])
Fig. 1. Terra satellite, launched 18th December 1999 (orbit height: 713 km; orbital velocity: 7 503 km·s-1; maximum velocity: 27 010 km·h-1) and the position of the MODIS sensor (source: [8])

Snímky označené identifikátorem MOD10_L2 poskytují informaci o sněhové pokrývce v denním kroku. Detekce probíhá pomocí normalizovaného diferenčního sněhového indexu (Normalized Difference Snow Index, NDSI). Dalším produktem je série korekčních snímků určených ke zmírnění chyb a označení detekce nejisté sněhové pokrývky. Zasněžená krajina má obvykle velmi vysokou odrazivost ve viditelných pásmech a velmi nízkou odrazivost pro krátkovlnná infračervená pásma. NDSI odhaluje velikost tohoto rozdílu. Každá datová granule obsahuje 5 minut dat ze svazku pozorovaných v rozlišení 500 m. Sběr dat byl zahájen 24. února 2000 a v současnosti probíhá revize dat pro aktuální verzi 61, která by měla být dokončena na jaře 2022.

Stažení družicových dat a jejich zpracování v prostředí GIS

Praktická stránka zpracování těchto dat v ČHMÚ spočívá v přednastavení parametrů zájmové oblasti na obdélník překrývající území Česka. Pokud je pořízený snímek v průniku s tímto obdélníkem, jsou v nejbližším termínu pracovníci ČHMÚ e-mailem informováni o jeho dostupnosti s odkazem ke stažení daného datového souboru. Termín zaslání závisí na komplexnosti snímku. V případě, že se snímek vyznačuje velkým množstvím tříd, může docházet ke zpoždění. Nejčastěji však bývají notifikace zasílány do 12 hodin od pořízení.

Data jsou poskytována ve formátu HDF-EOS2 a jsou ukládána jako 8bitová celá čísla bez znaménka. Datový formát HDF (Hierarchical Data Format) umožňuje efektivně ukládat rozsáhlá, a přitom poměrně rozmanitá data a metadata [9]. Pro zájmovou oblast zahrnující území Česka je velikost takových souborů přibližně 5–25 MB, což odráží plochu a prostorovou distribuci sněhu v krajině. Každý HDF soubor je složen z několika parametrů, z nichž je pro detekci sněhu zásadní výstup „NDSI_Snow_Cover“, který obsahuje příznaky rozdělené do devíti tříd uvedených v tab. 1.

Tab. 1. Třídy výstupu NDSI_Snow_Cover
Tab. 1. Classes of the NDSI_Snow_Cover output

Pro navazující práci v prostředí GIS je nutné extrahovat jednotlivé třídy nejprve do rastrové podoby a následně do polygonů, z nichž je využita zásadní hranice mezi oblastí se sněhem a bez sněhu. V první fázi extrahování je nutné použít nástroj HEG (HDF-EOS To GeoTIFF Conversion Tool), který je volně dostupný jako podpůrný software z portálu NASA [10]. Umožňuje dostatečně přesnou konverzi z HDF formátu do GeoTIFF formátu tak, aby při projekci UTM a zadání odpovídající zóny (pro Česko 33N, nebo 34N na východě) došlo k plnohodnotnému překryvu s českým digitálním modelem reliéfu (DMR; k němu viz dále). Nejlepším ověřením kvality překryvu jsou lokální odrazy ve třídě 237, tedy vnitrozemské vodní plochy, které lícují s hydrografickým podkladem v GIS (např. vodní nádrže Rozkoš nebo Nové Mlýny). Příklad nastavení nástroje HEG znázorňuje obr. 2.

Obr. 2. Nastavení parametrů pro konverzi v nástroji HEG
Fig. 2. Parameter settings for conversion in the HEG tool

Zpracování snímků v prostředí ArcGIS Desktop lze rozdělit do několika fází, přičemž použití jednotlivých nástrojů se přizpůsobuje možnostem licence ArcGIS Desktop v rámci ČHMÚ a do budoucna lze též předpokládat adaptaci na novější verze a produkty:

  1. Extrakce zájmových dat ze snímku DPZ. Snímek ve formátu GeoTIFF je importován s parametry nastavenými v nástroji HEG (projekce, rozsah). Nejčastěji je pro oblast Česka pořízen snímek v dopoledních hodinách, jenž obsahuje území přibližně od jižní Skandinávie po oblast Alp. V méně častých případech lze použít i snímky jen s částečným překryvem, což vychází z orbitální dráhy družice, která území Česka může zaznamenat z více přeletů. Snímky, jež mají průnik s přednastavenou obdélníkovou maskou pro Česko, jsou podrobeny ořezu, aby se snížila náročnost dílčích výpočtů. Získaný rastr je nutné nejprve reklasifikovat s ohledem na vymezení oblastí se sněhem a bez sněhu. Původní data do roku 2017 byla koncipována obecněji do tří základních tříd: sníh, bez sněhu a oblačnost. Po roce 2017 je již kvalita odrazu ze sněhové pokrývky rozdělována detailněji do tříd o hodnotách 0 až 100 (dle indexu NDSI), kde 0 je zaručená oblast bez sněhu a 100 je nejvyšší možná odrazivost od sněhové pokrývky. Pro účely hydrologie ČHMÚ jsou veškeré hodnoty pro odraz sněhu v rozmezí 1 až 100 považovány jako sníh, aby bylo možné získat co největší balík dat pro vyhodnocení nulové izochiony. Ve výsledku jsou tedy všechny hodnoty 1–100 překlasifikovány do hodnoty 50 a hodnoty 0 jsou brány jako území bez sněhu. Významným prvkem vycházejícím z aktuálních podmínek počasí je třída o hodnotě 201, což je území se sporným vyhodnocením (interní výpočet NDSI), a třída o hodnotě 250, což je oblačnost, která je nejčastějším limitujícím faktorem pro plnohodnotné vyhodnocení. Ostatní třídy lze považovat za doplňkové a je možné je z dalších operací vynechat. Reklasifikovaný rastr je ještě generalizován nástrojem Boundary Clean, aby byly potlačeny dílčí mikroregiony o velikostech jednotek pixelů, a poté jsou třídy barevně označeny, aby vynikl primární vizuální přehled o datové sadě. Pro další operace je nejprve nutné převést rastrovou vrstvu na polygonovou a tu následně na liniovou pomocí nástroje Feature to Line. Takto vznikne datová sada linií, jež mají specifický gridcode dle původního rastru, který obklopovaly. Konkrétně, na kontaktu dvou polygonů, kde jeden je o hodnotě 0 (bez sněhu) a druhý o hodnotě 50 (sníh), vzniknou dvě linie, jedna s gridcodem 0 a druhá s gridcodem 50. Hledaným prvkem je následně vrstva, která vznikne průnikem výběrů (nástroj Intersect) těchto dvou linií a obsahuje všechny viditelné hranice mezi oblastmi se sněhem a bez sněhu. Tuto nesouvislou linii ohraničující zaznamenanou sněhovou pokrývku lze označit za nulovou izochionu. Tato linie již poskytuje jistou prostorovou představu o pozici hranice sněhu v rámci území Česka. Pro bližší údaj o pozici je nutné najít přibližnou hodnotu nadmořské výšky, v níž se tato linie nachází. Jako nejjednodušší metoda byla zvolena extrakce pixelů pomocí masky, kde za masku je považována právě linie izochiony. Takto jsou extrahovány hodnoty nadmořské výšky z rastrového podkladu, kterým je DMR o rozlišení 25 m (viz dále).
  2. Prostorová analýza dat. Jako nejvhodnější dělení území Česka byly pro potřeby definice nulové izochiony zvoleny geomorfologické oblasti. Ty nejlépe odrážejí vlastnosti reliéfu, relativní i absolutní členitost (pohoří/nížiny) a rozdělení orientace svahů k jednotlivým světovým stranám (sever/jih, západ/východ), tedy faktory, u nichž je předpokládán zásadní vliv na akumulaci a tání sněhové pokrývky. V rámci Česka je takových zpracovávaných oblastí 27 (blíže viz tab. 2 a obr. 3). Pro každou z těchto oblastí je pomocí nástroje Zonal Statistics as Table vyhodnocen extrahovaný soubor pixelů přináležející dané geomorfologické oblasti a výstupem je statistika obsahující informace o počtu pixelů, jejich minimální a maximální hodnotě, rozsahu hodnot, sumě hodnot a především průměrné hodnotě, což je položka, se kterou je následně pracováno. Vektorová vrstva geomorfologických oblastí supluje funkci databáze poloh nulové izochiony, jelikož pro každý analyzovaný den je vytvořen sloupec hodnot, kde je každé geomorfologické oblasti přiřazena hodnota průměrné polohy, pokud v daný den byla zaznamenána. Pro propojení statistického výstupu s atributovou tabulkou vrstvy oblastí je použit nástroj Join Field.
  3. Vizuální interpretace dat. Vizualizaci hodnot lze provést pomocí labelingu hodnot pro jednotlivé oblasti, případně i s podbarveným kartogramem. Situace, kdy je definována průměrná hodnota pro každou nebo alespoň pro většinu geomorfologických oblastí, je během pozorování spíše výjimečná, neboť se často projevuje faktor oblačnosti a při absenci oblačnosti je sníh buď omezen jen na horské oblasti, nebo naopak pokrývá celé území Česka. Jak již bylo zmíněno výše, hranice sněhu vychází z různé intenzity indexu odrazivosti NDSI a každou z průměrných hodnot je potřeba podrobit kritické analýze (neboli validovat), zda jde o objektivní hodnotu a zda může reprezentovat podmínky v dané oblasti. Při posouzení se vychází nejen z interních dat o sněhové pokrývce pocházejících z pozorovací sítě ČHMÚ (automatické stanice, pozorovatelé, terénní měření), ale například i z výstupů webových kamer nebo historické korelace mezi oblastmi. Do výpočtu zásob vody ve sněhu vstupuje poloha nulové izochiony jako limitní hodnota pro prostorovou interpolaci parametrů sněhové pokrývky, kdy je pro každou z oblastí vygenerována virtuální síť nulových bodů, jež zamezuje interpolaci odhadovat nenulovou (kladnou) výšku sněhové pokrývky nebo vodní hodnotu sněhu pod touto pozicí. Pro tuto analýzu v prostředí GIS se používá nástroj ClidataGIS, který umožňuje import naměřených dat, jejich detailnější vizuální kontrolu i nastavení parametrů pro interpolaci. Finálním výstupem je mapa včetně doplňující tabulky zohledňující změny ve výskytu sněhu v předchozích týdnech, jež je vyvěšena na portálu ČHMÚ, a veřejnost se zde může seznámit s předpokládaným objemem vody v dílčích zájmových povodích (vodní nádrže, významné závěrové profily vodních toků). Aktuálně je takový výstup, vycházející z pondělních změřených hodnot, generován jednou týdně, a to v úterý. Do budoucna bude ovšem díky automatizované síti sněhoměrných stanic v kombinaci se satelitními snímky možné zpracovat podobné analýzy i častěji během týdne.
Tab. 2. Geomorfologické oblasti, pro které je v ČHMÚ v zimní sezoně určována průměrná nadmořská výška nulové izochiony, a jejich identifikátory (upraveno podle [11])
Tab. 2. Geomorphological regions for which, in the winter season, the CHMI determines the average altitude of the zero isochion, and their identifiers (adapted from [11])
Obr. 3. Geomorfologické oblasti (pro ID viz tab. 2) a jejich nadřazené subprovincie v Česku (upraveno podle [11])
Fig. 3. Geomorphological regions (for ID see tab. 2) and their parent subprovinces
in Czechia (adapted from [11])

Základní geografické podklady pro získání terénních vysvětlujících proměnných

Na prostorovou (ale i časovou) variabilitu nadmořské výšky nulové izochiony má vliv konfigurace terénu. Časová variabilita bude jistě více spjata s klimatickými podmínkami územních celků, pro které je studium prováděno. Vzhledem k tomu, že těmito územními celky byly z výše uvedených důvodů geomorfologické oblasti, bylo nutné získat vektorovou vrstvu s polygony představujícími tyto regiony. Ta byla stažena z Geoportálu Českého úřadu zeměměřického a katastrálního (ČÚZK), kde je součástí databáze Data200 (konkrétně vrstva Popis) [12]. Tato vrstva vychází z geomorfologického členění obsaženého v publikaci [11], jež ve skutečnosti odkazuje na 28 oblastí. V ČHMÚ se však tradičně uvažuje pouze 27 oblastí, neboť Záhorská nížina je sloučena s Jihomoravskou pánví (řádek s ID XA v tab. 2). V tomto smyslu byla také vrstva geomorfologických oblastí před dalšími analýzami upravena.

Digitální model reliéfu (DMR) v podobě rastru, jenž byl zdrojem informací o nadmořské výšce a jiných terénních parametrech v geomorfologických oblastech, vychází z Digitálního modelu území zpracovaného v měřítku 1 : 25 000 (tzv. DMÚ 25), který ČHMÚ zakoupil v roce 2001 od Vojenského geografického a hydrometeorologického úřadu (VGHMÚř) generála Josefa Churavého. Tento rastr se čtvercovými buňkami o straně 25 m vznikl z původních dat přímo v ČHMÚ, přičemž zásadním pro jeho vytvoření byl výškopis ve formě vrstevnic. DMR tak vznikl vhodnou interpolací tehdy ještě v systému S-42. Ale protože ČHMÚ postupem času přešel na systém UTM zone 33/34N, byl rastr reprojektován a převzorkován právě do tohoto systému (přesněji pouze pásma 33). Stejně tak byla i geometrie vrstvy polygonů geomorfologických oblastí transformována do systému UTM zone 33N.

Za využití polygonové vrstvy a rastru následovalo extrahování několika charakteristik terénu pro každou geomorfologickou oblast. Ty měly sloužit jako vysvětlující proměnné v plánované regresní analýze. Jejich výčet a význam uvádí tab. 3. Přitom je nutné speciálně upozornit na ukazatel SDASP, který v české literatuře není příliš znám. Jde o tzv. směrovou směrodatnou odchylku odrážející variabilitu orientace svahů v jednotlivých regionech vyjádřenou v radiánech [13]. Pro výpočet charakteristik terénu souvisejících s úhly byl využit R balíček circular [14]. S vektorovou vrstvou bylo při extrakci manipulováno prostřednictvím R balíčku sf [15], s rastrem naopak pomocí R balíčku terra [16].

Tab. 3. Charakteristiky geomorfologických oblastí Česka získané jako potenciální vysvětlující proměnné pro návaznou regresní analýzu
Tab. 3. Characteristics of geomorphological regions of Czechia obtained as explanatory variables for next regression analysis

Statistické zpracování získaných dat o nadmořské výšce nulové izochiony

Obr. 4 prozrazuje, kolik hodnot nadmořských výšek nulové izochiony bylo v původní databázi k dispozici od prosince 2012 do května 2021, kdy již pracovníci ČHMÚ vyhodnocovali data kompletně svépomocí. Informace o nadmořské výšce nulové izochiony byla přirozeně k dispozici pouze za zimní období, která obyčejně bývají rozdělena do dvou kalendářních let. Toto rozdělení však není pro další zpracování časových řad příliš vhodné, a proto byly zimní sezony nadále přiřazovány k rokům s větším podílem měsíců. Jelikož za zimní sezonu bylo považováno období listopad až květen, měsíce listopady a prosince byly přiřazeny k následujícím kalendářním rokům.

Obr. 4. Celkový počet dní (hodnot) s detekovanou nulovou izochionou za všechna období zim 2013–2021
Fig. 4. Total number of days (values) with the detected zero isochion for all winter
seasons 2013–2021

Ve výsledku tedy mohly být zimní sezony analyzovány za období 2013–2021, tedy za devět let. Z nejrůznějších důvodů popsaných již dříve nebyla původní databáze pro jednotlivé geomorfologické oblasti kompletní, přičemž dokonce nebyl dodržen ani ekvidistantní týdenní krok, neboť např. z důvodu oblačnosti musel být vybrán některý z následujících bezoblačných dnů v týdnu. Chybějící hodnoty a nedodržení stejného časového kroku činí taková data pro následující statistické zpracování poměrně problematickými, protože naprostá většina statistických modelů vyžaduje úplnost a konstantnost časového kroku (zejména jde-li o modely časových řad). Přestože se vyvíjejí také modely pro data tohoto typu (např. [17]), obecně je doporučováno zbavit se jejich výše zmíněných nedostatků, aby mohly být aplikovány tradiční modely. Databáze hodnot nadmořských výšek nulové izochiony byla nakonec doplněna odhady hodnot tak, že pro každou geomorfologickou oblast vznikla časová řada v denním kroku, jež bez jediné chybějící hodnoty reprezentovala všechny zimní sezony období 2013–2021. Tohoto výsledku bylo dosaženo prostřednictvím modifikovaného EM (Expectation-Maximization) algoritmu, který s neúplnou časovou řadou zachází jako s řadou vícerozměrnou, kde je tedy uvažován vektor jednorozměrných řad s vazbami mezi nimi a zároveň je na všechny jeho prvky jako filtr aplikována metoda splinů. Mnohem podrobněji je algoritmus popsán v [18], přičemž stejní autoři jej implementovali do R balíčku mtsdi [19], jehož funkce mnimput byla také využita pro doplnění chybějících denních hodnot.

Denní hodnoty byly nadále agregovány pro jednotlivé roky (respektive zimní sezony) pomocí alfa-useknutého průměru, aby se předešlo citlivosti na extrémy (blíže k jeho vlastnostem viz např. [20]). Useknuto bylo 10 % extrémních hodnot. Pro roky, u nichž to bylo možné, byly navíc hodnoty agregovány tak, aby reprezentovaly kromě kompletního období (značeného jako KOMPLET) také jednotlivé etapy zimních sezon (tj. AKUMULACE, TÁNÍ a nerozlišitelný ZBYTEK mezi akumulací a táním). Pro definici období akumulace, respektive tání, byl jako referenční použit vývoj sněhové pokrývky v pohraničních horách, především pak v Krkonoších a Jizerských horách. Nebylo přitom počítáno s možnými překryvy těchto období způsobenými případnými oblevami. Avšak různé trvání jednotlivých období v různých letech bylo zaručeno. Poté došlo taktéž k agregaci přes všechny roky, aby byly získány „dlouhodobé“ hodnoty jako vysvětlované proměnné pro připravované regresní modely. Datumy vymezující období akumulace a tání sněhu přirozeně nebyly ve všech letech stejné, a proto se postupovalo striktně podle toho, co bylo pozorováno kombinací družicových snímků a terénního průzkumu. Pro období TÁNÍ nebylo datum jeho začátku v případě několika geomorfologických oblastí k dispozici, a proto bylo nahrazeno posledním datumem období ZBYTEK. Překryvy období uvažovány nebyly, spíše bylo cíleno na delší charakter období.

Z několika důvodů, mezi něž patřilo i množství získaných dat determinované počtem geomorfologických oblastí, byl pro samotnou regresní analýzu nakonec vybrán lineární model (vícerozměrný aditivní s odhadem parametrů pomocí metody obyčejných nejmenších čtverců) a model náhodných lesů. Terénní vysvětlující proměnné byly v prvním případě vybírány na základě Akaikeho informačního kritéria (za kombinace dopředného a zpětného hledání proměnných; blíže viz [21, 22]). Výběr vysvětlujících proměnných konkrétně probíhal prostřednictvím funkce stepAIC implementované v R balíčku MASS, jenž je součástí knihy [23]. Před aplikací lineárních modelů byla speciálně zkoumána možná kolinearita mezi nabízenými vysvětlujícími proměnnými pomocí Pearsonových korelačních koeficientů. V případě náhodných lesů byl aplikován dopředný výběr pomocí funkce ffs implementované v R balíčku CAST [24–26], který ke svému fungování potřebuje především R balíčky caret [27, 28] a randomForest [29].

Ke zjištění, zda lze v nadmořské výšce nulové izochiony pozorovat významný meziroční monotónní trend, byl pro jednotlivé geomorfologické oblasti aplikován neparametrický Mannův–Kendallův test. Vzhledem k tomu, že tento test je i přes svůj neparametrický charakter citlivý na autokorelovanost v časové řadě, byla pro řadu období KOMPLET (2013–2021) zvolena modifikace TFPW (trend-free pre-whitening) implementovaná v R balíčku zyp [30]. Teoretické základy této modifikace testu lze studovat v publikacích [31–33]. Ostatní části zimních období nemohly být takto studovány, protože získané časové řady jsou zatím velmi krátké.

VÝSLEDKY A DISKUZE

Obr. 4. ukazuje, kolik hodnot nadmořských výšek nulové izochiony se celkem podařilo ze snímků MODIS odvodit pro každou z 27 geomorfologických oblastí používaných v praxi hydrologie ČHMÚ. Celkově bylo tedy za období prosinec 2012 až květen 2021 k dispozici 3 216 takových hodnot, což z možného celkového množství denních hodnot (při uvážení všech zimních období obsažených v době od listopadu 2012 do května 2021, tj. 51 570 hodnot) činí zhruba 6,24 %. Není to způsobeno jen tím, že se hydrologové ČHMÚ tradičně soustředí pouze na získání jedné hodnoty týdně pro každou oblast, ale také tím, že některé oblasti postihuje oblačnost mnohem častěji než jiné. Zároveň je třeba podotknout, že u nížinných oblastí je pravděpodobnost zasněžení mnohem menší než u oblastí vyznačujících se spíše horským reliéfem. Na obr. 4 je tento fakt zdůrazněný metodou nepravého kartogramu, kde jsou geomorfologické oblasti zařazeny do jednotlivých tříd, které se vyznačují různou intenzitou modré barvy. Rozdíl, jenž je zcela jistě způsoben typickou nadmořskou výškou nebo členitostí uvnitř regionů, je velmi dobře patrný. Kromě výše uvedeného je nutné zohlednit i fakt, že vzhledem k plošnému výskytu i době setrvání sněhové pokrývky lze zimy v předešlých pěti letech hodnotit spíše jako podprůměrné, což lze názorně sledovat na průběžných statistikách [2].

V tab. 4 jsou výsledky týkající se extrakce vybraných terénních proměnných, u kterých bylo předpokládáno, že mohou vysvětlit prostorovou variabilitu nadmořské výšky nulové izochiony, a mohou být tedy významné pro sestavení regresních modelů. Seznam proměnných, jež lze takto odvodit z DMR, není jistě vyčerpávající, ale předpokládalo se, že reprezentuje alespoň ty nejdůležitější faktory související se zeměpisnou šířkou, délkou a střední, minimální i maximální nadmořskou výškou, členitostí a sklonitostí terénu, stejně jako i s orientací svahů. Zvláštní postavení zde zaujímají charakteristiky týkající se variability nadmořské výšky, ale též orientace svahů vůči světovým stranám. Ukazatele spojené s orientací svahů bylo možné vyjádřit v úhlech, ale nakonec bylo rozhodnuto, že do regresních modelů vstoupí jako kategorické proměnné vzniklé reklasifikací, neboť v takových modelech není jednoduché počítat s úhlovými veličinami a tak jako tak je doporučována jejich transformace, přičemž jsou zpravidla aplikovány goniometrické funkce.

Tab. 4. Hodnoty vybraných terénních vysvětlujících proměnných pro 27 geomorfologických oblastí Česka
Tab. 4. Values of selected explanatory variables related to terrain of the 27 geomorphological regions of Czechia

Tab. 5 poskytuje alespoň základní představu o tom, jaké jsou dlouhodobé (lépe řečeno dlouhodobější) hodnoty nadmořské výšky nulové izochiony na území jednotlivých geomorfologických oblastí Česka. Kromě toho je v tab. 5 ukázáno, jak se tyto dlouhodobé charakteristiky liší podle jednotlivých etap zimní sezony (tj. období akumulace, tání a prostřední období, kdy akumulaci či tání nelze rozlišit), a to alespoň pro roky 2018–2021. Je třeba poznamenat, že úplná vícerozměrná řada mohla díky extrapolacím obsahovat též hodnoty větší než maximální a menší než minimální nadmořské výšky vyskytující se v geomorfologických oblastech, což v doplněných datech nastalo u 16 oblastí (jinými slovy u 1,3 % z celkového počtu denních hodnot). To ale při následujícím použití modelů pro zjišťování, které terénní charakteristiky mají vliv na variabilitu nadmořské výšky nulové izochiony, překážku nekladlo. Navíc tyto hodnoty by bylo možné považovat za reálné, kdyby např. nejvyšší partie oblastí dosahovaly odhadovaných pozic. Je rovněž předpokládáno, že tyto situace byly redukovány aplikováním alfa-useknutého průměru.

Tab. 5. Dlouhodobé průměry nadmořské výšky nulové izochiony za zimní období (listopad–květen) pro 27 geomorfologických oblastí Česka
Tab. 5. Long-term averages of the zero isochion altitude for winter seasons (November–May) for the 27 geomorphological regions of Czechia

Před aplikováním lineárních modelů je doporučováno všímat si zejména pravděpodobné kolinearity mezi vysvětlujícími proměnnými, tj. jevu, při němž dvě nebo více proměnných poskytují velmi podobnou informaci. Obr. 5 prostřednictvím Pearsonových korelací ukazuje, že i přes neúplný seznam terénních proměnných velmi pravděpodobně kolinearita v sadě vysvětlujících proměnných přítomna byla. Zejména si lze povšimnout velmi těsné (a statisticky významné korelace na hladině 0,05) mezi maximem nadmořské výšky a rozpětím mezi minimem a maximem. Dále je možné vypozorovat velmi těsný vztah mezi oběma veličinami vztahujícími se k orientaci svahů. Z těchto důvodů nebylo v případě lineárních modelů dále počítáno s rozpětím mezi maximem a minimem nadmořské výšky a s převažující orientací svahů. Alternativou (při ponechání všech dostupných vysvětlujících proměnných) mohly být regresní modely, v nichž jako vysvětlující proměnné figurují namísto původních proměnných hlavní komponenty (viz např. [34]). U modelů náhodných lesů byly před výběrem záměrně ponechány všechny získané terénní vysvětlující proměnné. Poznamenejme ještě, že pro modely náhodných lesů bylo ponecháno původní nastavení parametrů jako např. v [35].

Obr. 5. Pearsonovy korelace mezi vybranými terénními vysvětlujícími proměnnými (elipsami a barvami znázorněny pouze korelace významné na hladině 0,05)
Fig. 5. Pearson’s correlations between selected explanatory variables related to terrain (by ellipses and colours depicted only correlations significant at the 0.05 level)

Tab. 6 a 7 již prozrazují, které terénní proměnné byly konkrétně vybrány pro lineární modely (za asistence Akaikeho informačního kritéria), respektive pro modely náhodných lesů (za pomoci dopředného algoritmu podle [25]). Je vidět, že lineární modely jsou mnohem komplexnější, pokud jde o inkluzi vysvětlujících proměnných. Některé proměnné dle statistiky t významné nejsou, ale i tak přispívají k významnosti celých modelů. Ty jsou dle statistiky F významné dokonce na menších hladinách než 0,05. Rovněž tak hodnoty koeficientů determinace (R2) jasně naznačují, že vysvětlení variability dlouhodobé nadmořské výšky nulové izochiony je zde více než dobré. Za cenu snížení hodnot R2 algoritmus u modelů náhodných lesů vybral vysvětlujících proměnných méně. Velmi často zde vystupuje zeměpisná délka a extrémy nadmořské výšky. Stálou vysvětlující proměnnou je tu charakteristika související se střední nadmořskou výškou geomorfologických oblastí, což potvrzuje situaci na obr. 4. Pro období akumulace a tání sněhu se zdají být důležitými také charakteristiky spjaté s varia-
bilitou nadmořské výšky (v případě akumulace) a orientace svahů (v případě tání), což zní poměrně logicky. Je však třeba poznamenat, že náhodný les je model založený na resamplovacích technikách, takže při jiném běhu algoritmu může dojít k nepatrně odlišnému výběru proměnných. Domníváme se ale, že i tak by tyto výběry byly velice podobné. Např. pro tání sněhu bude důležitá variabilita orientace svahu, jež souvisí s příznivými nebo naopak nepříznivými podmínkami v průběhu světlé části dne.

Tab. 6. Nejlepší lineární modely podle Akaikeho informačního kritéria
Tab. 6. Best linear models according to the Akaike Information Criterion

t – kvantil Studentova t-rozdělení / Student t-distribution quantile; P – pravděpodobnost / probability; F – kvantil Fisherova-Snedecorova F-rozdělení / Fisher-Snedecor F-distribution quantile; p – p-hodnota / p-value

Tab. 7. Vysvětlující proměnné vybrané algoritmem podle [25] pro modely náhodných lesů
Tab. 7. Explanatory variables selected by the algorithm of [25] for random forest models

Analýza trendů, a tedy meziroční časové variability nadmořské výšky nulové izochiony byla provedena jen pro nejdelší časovou řadu označenou jako KOMPLET, protože období čtyř let, pro které je k dispozici rozdělení podle jednotlivých fází přírůstku a úbytku sněhové pokrývky, není možné pro tento typ analýzy ještě považovat za reprezentativní. Z tab. 8 je jasné, že během zimních sezon 2013–2021 byla nulová izochiona spíše stabilní. Přesto je však možné si povšimnout, že v pěti regionech pravděpodobně dochází k poklesu (Podkrušnohorská oblast, Středočeská tabule, Severní Vněkarpatské sníženiny, Jihomoravská pánev) či vzestupu nulové izochiony (Krkonošsko-jesenické podhůří). Důvody pro vznik takových trendů mohou být různé, od skutečných nárůstů či poklesů sněhové pokrývky až po fakt, že data pro některé z těchto regionů nemusela být dostatečná.

Tab. 8. Výsledky trendové analýzy pro všech 27 geomorfologických oblastí Česka ( : statisticky významný rostoucí trend na hladině 0,01; : statisticky významný klesající trend na hladině 0,01)
Tab. 8. Results of trend analysis for all 27 geomorphological regions of Czechia ( : statistically significant increasing trend at the 0.01 level; : statistically significant decreasing trend at the 0.01 level)

Tak např. při pohledu na obr. 6, kde jsou černou čarou znázorněny průběhy časových řad jen pro regiony se signifikantním monotónním trendem, je zřejmé, že nejméně tři z těchto výsledků jsou dosti nevěrohodné. K situacím, kdy R2 = 1, nedochází téměř nikdy. Navíc průběhy těchto podezřelých řad nevykazují téměř žádnou variabilitu (např. černé čáry jsou zakryty modrými regresními přímkami), což naznačuje, že spíše selhal EM algoritmus při doplňování chybějících hodnot, který pracoval pouze s několika málo zjištěnými hodnotami, jež mohly být nadto zatíženy velkou nejistotou. Doplňme, že zatímco v obr. 6 jsou regresní přímky konstruovány za využití lineárních modelů, tak v tab. 8 jsou regresní koeficient a absolutní člen pro porovnání výsledků vztaženy k tzv. Senovu neparametrickému odhadu [36].

Obr. 6. Průběh ročních řad průměrné nadmořské výšky nulové izochiony (2013–2021) v pěti geomorfologických oblastech Česka, u kterých byl nalezen statisticky významný monotónní trend na hladině 0,01 (černá čára: konkrétní hodnoty výšky; modrá přímka: lineární trend získaný metodou obyčejných nejmenších čtverců)
Fig. 6. Course of the annual series of average zero isochion altitude (2013–2021) in five geomorphological units in Czechia for which a statistically significant monotonic trend was found at the 0.01 level (black line: specific altitude values; blue line: linear trend obtained by the ordinary least squares method)

Pro účely aktualizace tabulky diferencí mezi nadmořskými výškami nulové izochiony v jednotlivých geomorfologických oblastech, kterou uvádí [3], byla sestavena její nová verze (viz tab. 9). Je patrné, že současná čísla se dosti odlišují od těch publikovaných v minulosti. Mohlo např. dojít k určitému zpřesnění, kde aplikovaný alfa-useknutý průměr již k extrémním hodnotám přihlížel jinak. Z praktického hlediska metodiky výpočtu zásob vody ve sněhu je však pro prognostiky ČHMÚ nejzásadnější vztah mezi horskými pohraničními regiony, kde se sněhová pokrývka vyskytuje nejdéle a nejčastěji. Nutné je ale zmínit i riziko významných jarních povodní z tajícího sněhu, jež je spjato spíše se sněhovou pokrývkou v nížinách, tzn. v tabulích v Polabí. V takových situacích je nulová izochiona buď zcela potlačena a sníh se vyskytuje na celém území, nebo je hranice poměrně ostrá a ohraničuje nejteplejší oblasti a tepelné ostrovy. Tab. 9 zkušenosti z analýz nulové izochiony jen potvrzuje, především pak nejpoužívanější vztah, a to ten mezi Krkonošskou oblastí a Šumavskou hornatinou. Ten je možné charakterizovat zjednodušeně tím, že sněhová pokrývka na Šumavě začíná vždy o 100 a více výškových metrů výše oproti hranici sněhu v Krkonoších. Důvody pro to lze spatřovat v členitosti místních regionů, kdy Krkonoše reprezentují spíše strmé svahy a českou část Šumavy naopak pozvolnější přechod do nižších oblastí.

Tab. 9. Rozdíly mezi dlouhodobými průměry nadmořské výšky nulové izochiony za zimní období (listopad–květen) napříč všemi 27 geomorfologickými oblastmi Česka
Tab. 9. Differences between long-term averages of the zero isochion altitude for winter seasons (November–May) across all 27 geomorphological regions of Czechia

Obr. 7 demonstruje příklad výstupu získaného pomocí snímků MODIS a zároveň dokumentuje podmínky úbytku sněhové pokrývky začátkem dubna 2021. Je zde zobrazena nejčastější situace, kdy je sněhová pokrývka přítomna jen v nejvyšších polohách hor, kde je navíc částečně překryta oblačností. Naopak v nížinách už sníh neleží vůbec, stejně jako se na konci zimy již nevyskytují nejistě definovatelná území se sněhem, jež charakterizuje kód 201 (viz tab. 1).

Obr. 7. Analyzovaná situace z 9. dubna 2021, kde jsou vybrané třídy dle tab. 1 zobrazeny nad podkladem DMR (hodnoty určují průměrnou polohu nulové izochiony v m n. m. a 0 znamená oblast bez zaznamenané nulové izochiony)
Fig. 7. Analyzed situation on 9th April 2021 where the selected classes from tab. 1 are depicted above the digital elevation model (values determine the average position of the zero isochion in m a.s.l. and 0 determines the area without recorded zero isochion)

ZÁVĚR

V příspěvku jsou uvedeny aktualizované poznatky týkající se nadmořské výšky nulové izochiony, tedy čáry představující hranici mezi prostorem se sněhem a prostorem bez sněhu, ve 27 geomorfologických oblastech Česka, které se v praxi hydrologie ČHMÚ využívají pro odhad vodní hodnoty sněhu. Jelikož sníh představuje významnou komponentu tvořící odtok na území Česka, je tato aktivita nezbytná před spouštěním hydrologických modelů používaných v ČHMÚ jak pro operativní účely, tak i pro bilanční výpočty. Terénní průzkum je finančně i časově nákladný, a navíc neposkytuje dostatečně detailní představu o prostorovém rozložení sněhu v regionech, v Česku často výškově dosti členitých. Proto se pro zpřesnění této představy úspěšně využívají družicové snímky, které po vhodné kalibraci pomocí dat sebraných v terénu mohou být velmi nápomocny i při určení hranice mezi regionem se sněhem a bez sněhu (viz např. [1] a [37]). Studie tak navázala na předchozí výzkum a za využití prodloužené časové řady nadmořských výšek nulové izochiony do května 2021 se snažila odpovědět na otázky nastíněné v závěrech práce [3]. Především byla studována prostorová variabilita nadmořské výšky nulové izochiony a prostřednictvím regresních modelů bylo zjišťováno, na kterých faktorech souvisejících s terénem (určených v prostředí GIS z DMR vycházejícího z díla DMÚ 25) závisí kolísání nulové izochiony nejvíce. Za reprezentativní vysvětlovanou proměnnou zde byl brán alfa-useknutý průměr vypočtený z denních hodnot, u něhož bylo předpokládáno, že redukuje vliv nejistoty v určení polohy nulové izochiony. Průměr nadmořské výšky nulové izochiony v jednotlivých regionech navíc do určité míry představuje limitní polohu, pod kterou již nelze v interpolačních procesech očekávat nenulovou (kladnou) hodnotu výšky sněhové pokrývky a potažmo i vodní hodnoty sněhu. Vzhledem k nejistotám spjatým se snímky MODIS se domníváme, že pokud bude odvozování nulové izochiony vhodně kombinováno s interpolací hodnot získaných pozemním měřením, může část těchto nejistot být do značné míry redukována. Zde sestavené regresní modely vysvětlující prostorovou variabilitu mohou být tedy velmi užitečné např. při odhadu (průměrné) pozice nulové izochiony v situacích, kdy jsou regiony z převážné většiny zakryté oblačností. Bylo zjištěno, že velký vliv mají střední, ale i extrémní hodnoty nadmořských výšek. Svůj podíl má zřejmě i zeměpisná délka. V případě akumulace sněhu se k vysvětlujícím proměnným v regresních modelech přidává také směrodatná odchylka vypočítaná z nadmořských výšek vyskytujících se v dané geomorfologické oblasti. V případě tání sněhu se zdá být významná směrodatná odchylka získaná ze všech úhlů určujících orientaci svahů. Studium časové variability nadmořské výšky nulové izochiony bylo provedeno trendovou analýzou s cílem zjistit, zda lze vypozorovat přítomnost monotónních deterministických trendů. Aby byl redukován vliv autokorelace, byla přitom aplikována TFPW modifikace Mannova–Kendallova testu. Z výsledků vyplývá, že v průběhu zimních sezon se sněhová pokrývka, a tím i nulová izochiona, chovala stabilně. Pouze v pěti geomorfologických jednotkách byl nalezen statisticky významný meziroční trend. Výsledky analýz je ale třeba interpretovat opatrně. Tak např. u tří oblastí vyšel koeficient determinace podezřele roven jedné, což je jev velmi vzácný. Spíše než na skutečnost tento fakt poukazuje na to, že při doplňování chybějících hodnot nadmořské výšky nulové izochiony využitý EM (Expectation-Maximization) algoritmus v některých případech selhal. Důvodem byl nedostatek informací, což bylo jistě způsobeno malým množstvím dat získaných ze snímků MODIS, respektive odvozených produktů.

Jak bylo naznačeno výše, na jaře 2022 má dojít k dokončení revizí produktu vycházejícího z indexu NDSI. Tím by mohlo v ČHMÚ následně proběhnout jakési ověření, zda jsou vztahy představené v tomto příspěvku stále platné. Přirozeně se nabízí nasazení (polo)automatického zpracování dat pocházejících ze snímků MODIS za předpokladu, že pracovníci ČHMÚ prohloubí své skriptovací znalosti. Velmi doporučováno je prodloužení časové řady směrem do historie s uvážením, že produkty spojené s detekcí sněhu na zemském povrchu jsou k dispozici již od roku 2000. Nabízí se také „zahuštění“ dat s ohledem na skutečnost, že snímky MODIS mají mnohem jemnější časový krok než jeden týden. To může umožnit i aplikaci sofistikovanějších modelů časových řad, než je pouhá analýza trendu. Zpřesní se tak představa o časové dynamice nulové izochiony, pokud i na straně vysvětlujících proměnných nebudou jen charakteristiky terénu, ale i řady klimatologických prvků. Vysvětlení rozdílů v poloze nulové izochiony lze totiž hledat též v rozdílu klimatických podmínek, především ve výskytu srážkových situací, jež ovlivňují převážně jihozápad území Česka, přičemž směrem k severním pohořím již nedosahují takových intenzit. Za úvahu rovněž stojí podívat se na další charakteristiky související s nulovou izochionou, které mohou do modelů vstupovat jako vysvětlované proměnné. Zajímavou informaci může hydrologům poskytnout např. gradient změny nulové izochiony nebo změna délky setrvání nulové izochiony nad určitou výškovou hranicí. Podobné ukazatele byly studovány již v [38] a zejména na období tání sněhu a jeho vliv na odtok byla zaměřena studie [39]. Za zmínku stojí rovněž metoda, jíž je v těchto publikacích nulová izochiona extrahována. Tato metoda se totiž liší od metodiky prezentované zde, a proto se tu přirozeně nabízí srovnání obou přístupů za využití dat pro území Česka.

Zásadním limitem pro přesnější definici polohy nulové izochiony zůstává prostorové rozlišení analyzovaného gridu (500 m). Pro hlubší implementaci nulové izochiony do výpočtu zásob vody ve sněhu je tedy vedle automatizace procesů nutné i prostorové zpřesnění. Jednou z možností by mohla být diferenciace reklasifikace produktů s hodnotami indexu NDSI podle různých faktorů, která se dle příspěvku [40] osvědčila na území Rakouska. Podmínkou však je, že délka odvození nulové izochiony bude vyhovovat požadavku operativní hydrologie ČHMÚ mít tento výsledek do 24 hodin. Aktuálně existují i přesnější prostorová data, ovšem opět se značným časovým prodlením mezi jejich pořízením a dostupností. Jistý příslib v horizontu této dekády lze spatřovat v lokálních pozorovacích systémech či v evropském programu Copernicus s misemi Sentinel (viz např. [41]). Neméně důležitým úkolem bude zajisté studium vlastností sněhové pokrývky separátně pro fáze akumulace, tání a ve zbytku zimní sezony, protože – jak výsledky této studie ukázaly – toto dělení má smysl.

V neposlední řadě se naskýtá otázka, zda není prostřednictvím snímků MODIS výhodnější mapovat plochu sněhové pokrývky přímo. To by bylo možné pouze při ideálních podmínkách za bezoblačnosti, kdy je zaručen přehled o celém území Česka, nebo alespoň o všech oblastech se sněhovou pokrývkou. Takových situací je ve skutečnosti pouze minimum, řádově jednotky dní během sezony. Z tohoto důvodu nemůže být aktuálně snímkování MODIS implementováno do výpočtů plochy se sněhovou pokrývkou jako pravidelně používaný nástroj, ale jen jako doplněk k informacím získaným z poměrně husté sněhoměrné sítě ČHMÚ. Naopak jistý potenciál lze najít právě v kalibraci modelu vycházejícího z dostupných řad nadmořské výšky nulové izochiony a jeho schopnosti predikovat pro „neviditelná“ místa.

Poděkování

Oba autoři jsou podpořeni v rámci tzv. Dlouhodobé koncepce rozvoje výzkumné organizace (DKRVO) ČHMÚ. Práce O. Ledvinky statistického charakteru je navíc účelově podpořena Technologickou agenturou České republiky (projekt SS01020366 „Využití dat dálkového průzkumu Země pro posouzení negativních dopadů přívalových srážek“). Autoři za tuto podporu vyslovují své díky.

 

Příspěvek prošel lektorským řízením.