Analýza dat pomocí R a ArcGIS Pro - 1. část

31. 5. 2018
  • Desktopový GIS
  • Další
  • ArcGIS Pro 2.1

Tento článek navazuje na předchozí články, ve kterých jsme probrali stažení a instalaci RStudia, propojení ArcGIS Pro a R, první analytické úlohyobohacení dat. Dnes se pustíme do podrobnější analýzy dat, která jsme si během tohoto procesu připravili.

Nejprve si ukážeme výpočet kriminality s ohledem na počet obyvatel, což nám poskytne lepší porovnání mezi plochami s různou hustotou zalidnění. Pro takto jednoduchou úlohu stačí samozřejmě použít Kalkulátor polí v ArcGIS, nicméně pro statisticky robustní analýzu využijeme možností R. Použijeme metodu Empirical Bayes smoothing, která vezme počet obyvatel v každém šestiúhelníku jako míru spolehlivosti dat. To nám umožní přizpůsobit výsledky v oblastech s vyšší a nižší důvěrou.

Připojení dat do R

Díky R-ArcGIS bridge jsou data z našeho projektu dostupná i v RStudiu. Nejprve otevřeme náš projekt v ArcGIS Pro a poté spustíme RStudio. Dále je potřeba načíst balíček arcgisbinding a následně ještě ověřit, zda přemostění funguje korektně a zda R správně rozezná verzi ArcGIS Pro. Všechny příkazy v konzoli RStudia spouštíme stisknutím tlačítka Enter.

  • Do konzole RStudia zadáme - library(arcgisbinding) a následně arc.check_product().
  • Nyní můžeme přistoupit k otevření dat v R – příkazem arc.check_product(<plná cesta k souboru s obohacenými daty San_Francisco_Crimes_Enrich_Subset.shp>) a vše uložíme do proměnné enrich_df. Tedy celý příkaz bude vypadat následovně:
    enrich_df <- arc.open(path = 'I:/Tipy a triky/2018/data/San_Francisco_Crimes_Enrich_Subset.shp')

Tímto příkazem jsme přiřadili proměnné enrich_df objekt typu arc.dataset. Objekt obsahuje jak prostorová, tak atributová data a je možné s nimi dále pracovat.

Převod datové sady

Práci s objektem si ukážeme při vytvoření další proměnné, která bude obsahovat pouze vybrané atributy původní proměnné enrich_df. Využijeme k tomu funkci arc.select(), jejímiž parametry jsou vstupní proměnná a seznam vybraných atributů. Vše uložíme do proměnné enrich_select_df. Celý příkaz bude vypadat takto:

enrich_select_df <- arc.select(object = enrich_df, fields = c('FID', 'SUM_VALUE', 'TOTPOP10', 'MEDHINC_CY', 'RENTER_CY', 'MEDVAL_CY', 'N13_BUS', 'N37_BUS', 'NLCDfrstPt')).

Tento nový objekt obsahuje prostorová data, ale pouze s devíti atributy. S typem proměnné data.frame dokáže R a všechny další komponenty pracovat. Protože naše data jsou prostorová, převedeme tuto proměnnou na tzv. spatial data frame, tedy na reprezentaci jednoho z následujících typů: SpatialPointsDataFrame, SpatialLinesDataFrame nebo SpatialPolygonsDataFrame (podle typu vstupních dat). Než však budeme pokračovat, je třeba do R nainstalovat balíček sp: install.packages("sp") a následně jej načíst do konzole příkazem: library(sp).

Vlastní převod proměnné provedeme příkazem
enrich_spdf <- arc.data2sp(enrich_select_df)
a dále tak budeme pracovat s proměnnou enrich_spdf.

Tímto krokem jsme ukončili napojení dat do R tím, že jsme je načetli a převedli do formy srozumitelné pro R.

Kam kliknout

Pokud chcete vyvolat nápovědu, v pravém sloupci je možné přepnout spodní okno na Help nebo jakmile se při psaní otevře bublina s nápovědou, stačí stisknout klávesu F1.

V pravém sloupci pak můžeme vidět přehled proměnných. Pokud na některou z nich klikneme, v dalším okně se nám zobrazí její obsah.

Analýza v R

Jak je vidět, názvy atributů úplně nekorespondují s podstatou dat, která obsahují, takže před samotnou analýzou ještě přejmenujeme některé sloupce.

Nejdříve si vytvoříme vektor řetězců s názvem col_names. Každá hodnota představuje nové názvy sloupců:
col_names <- c("OBJECTID", "Crime_Counts", "Population", "Med_HomeIncome", "Renter_Count", "Med_HomeValue", "Grocery", "Restaurant", "Pct_Forest").

Vlastní přejmenování pak provedeme funkcí colnames(), kdy proměnné enrich_spdf přiřadíme vektor nových názvů:
colnames(enrich_spdf@data) <- col_names.

Nyní můžeme dvojklikem na jméno proměnné v pravém okně zkontrolovat, že se přejmenování provedlo korektně. Stejného náhledu, ale tentokrát v konzoli, dosáhneme také příkazem
head(enrich_spdf@data).

Pro vlastní výpočet pomocí Empirical Bayes smoothing použijeme funkci EBest(), která je obsažena v balíčku spdep. Nemáme-li tento balíček, tak jej opět pomocí příkazu install.packages("spdep") nejprve nainstalujeme a poté příkazem library (spdep) jej načteme.

Pokud se při načítání balíčku objeví chyba, místo kopírování ze schránky vepište kód ručně.

Funkce samotná má dvě vstupní proměnné: xn. Za n zadáme vektor počtu kriminálních činů a za x počet obyvatel. Pro zpřehlednění můžeme celý výpočet zapsat takto:

n <- enrich_spdf@data$Crime_Counts

x <- enrich_spdf@data$Population

EB <- EBest (n, x)

Z výstupů funkce nás zajímá část, která vrací vypočtenou míru kriminality – a tu uložíme do proměnné p. Tedy: p <- EB$raw

Máme tedy nyní vypočteny hrubé odhady. Ty je ale dobré „vyhladit“ na základě důvěry v data, která se odvíjí od počtu obyvatel v každém šestiúhelníku. Z proměnné EB získáme informace o rozdělení (konkrétně průměrrozptyl) a tyto hodnoty uložíme do proměnných a, b.

b <- attr(EB, "parameters")$b

a <- attr(EB, "parameters")$a

Výsledný parametr, se kterým budeme dále pracovat, je proměnná z. Tu vypočteme tak, že odečteme před chvíli vypočtený průměr b od každé hrubé hodnoty p a výsledek vydělíme směrodatnou odchylkou.

Výsledný kód ve správné posloupnosti příkazů je tento:

b <- attr(EB, "parameters")$b

a <- attr(EB, "parameters")$a

v <- a + (b/x)

z <- (p - b)/sqrt(v)

Nakonec tato data přidáme jako nový atribut do objektu typu spatial data frame: vytvoříme nový atribut s názvem EB_Rate  a přiřadíme do něj hodnoty proměnné z.

enrich_spdf@data$EB_Rate <- z

Převod dat pro ArcGIS

Tím jsme zatím skončili s prací v R a vrátíme se do ArcGIS Pro. Objekt typu spatial data frame proto převedeme na objekt typu data frame.

arcgis_df <- arc.sp2data(enrich_spdf)

Tento objekt pak můžeme převést do SHP, tabulky nebo třídy prvků v geodatabázi pomocí funkce arc.write(). V našem případě zapíšeme data jako třídu prvků San_Francisco_Crime_Rates do souborové geodatabáze projektu ArcGIS Pro.

arc.write('C:/<cesta>/Analýza kriminality.gdb/San_Francisco_Crime_Rates', arcgis_df, shape_info = arc.shapeinfo(enrich_df)).

Vizualizace v ArcGIS Pro

V Katalogovém okně v ArcGIS Pro klikneme pravým tlačítkem na databázi projektu a zvolíme Obnovit. Objeví se třída prvků, kterou zapsalo R, a přidáme ji do nové mapy. Nyní provedeme hot spot analýzu s těmito daty, která by nám mohla ukázat oblasti s vysokou a nízkou mírou kriminality.

  • Na kartě Geoprocessing vyhledáme nástroj Optimalizovaná hot spot analýza a spustíme jej. Vstupní třídou budou San_Francisco_Crime_Rates a výstupní třídu nazveme jako San_Francisco_Crime_Rates_OHSA.shp. Polem pro analýzu bude atribut EB_Rate.

Jasně červené šestiúhelníky znázorňují oblasti, kde i po započtení vlivu počtu obyvatel je neobvykle vysoká míra kriminality. Zaznamenáníhodné je i to, že zde nyní nenajdeme oblasti s nízkou mírou kriminality (cold spots).

V dalším článku naši analýzu dokončíme tvorbou korelační matice.