Latvijas novadu kartogrammas veidošana ar paketi ggplot2

ggplot2 pakete programmā R ir ne tikai lieliska, lai veidotu dažādas diagrammas, bet izmantojot šo sistēmu ir iespējams veidot arī dažādas kartes, tai skaitā kartogrammas (choropleth map).

Lai izveidotu šādu kartogrammu ar Latvijas novadu robežām un informāciju par tiem, ir nepieciešamas divas lietas: (1) novadu robežas; (2) dati, ko vēlamies attēlot.

Novadu robežu iegūšanai var izmantot SIA Envirotech publiskoto GIS Latvija 10.2 datu bāzi, kas ir ESRI File geodatabase 10.2 datu formātā. Šo datubāzi var lejupielādēt no SIA Envirotech mājaslapas. Pēc lejupielādas nepieciešams atarhivēt zip failu.

Programmā R informāciju no šīs datubāzes var ielasīt ar speciālām funkcijām, kas pieejamas paketē rgdal. Pielāgojot piemēru, kas atrodas http://gis.stackexchange.com mājaslapā Aaron atbildē, darbības būs sekojošās:

Nārādām ceļu uz failu GIS_Latcija10.2.gdb, kas atrodas uz Jūsu datora.

library(rgdal)
Latvijagdb<-"/Users/didzis/Dropbox/DidzaDati/eksper/GIS_Latvija10.2/GIS_Latvija10.2.gdb"

Izmantojot funkciju ogrListLayers() iespējams apskatīties kādi ir slāņi datu bāzē.

ogrListLayers(Latvijagdb)

##  [1] "izobatas_line"            "izobatas_poly"           
##  [3] "izolinijas_line"          "krasta_linija_line"      
##  [5] "pagasti2001_poly"         "pagasti1993_poly"        
##  [7] "pagasti1999_poly"         "ramis_poly"              
##  [9] "robezas_poly"             "salas_poly"              
## [11] "udensteces_line"          "udenstilpes_poly"        
## [13] "pagasti2004_line"         "rajoni_line"             
## [15] "robezas_line"             "izolinijas_poly"         
## [17] "jurasNos_A0"              "jurasNos_A3"             
## [19] "gazes_vadi_line"          "ramis_line"              
## [21] "udenstilpes_line"         "apbuve_poly"             
## [23] "apbuve_point"             "ramis_mala"              
## [25] "purvi_poly"               "pagasti2004_poly"        
## [27] "pilsetas_point"           "rajoni_poly"             
## [29] "ciemi_point"              "dzelzcels_line"          
## [31] "novadi_2009_lin"          "novadi_2009_pol"         
## [33] "mezi_poly"                "jurasNos_A4"             
## [35] "Planosanas_reg_2009"      "pagasti2011_poly"        
## [37] "tks93_25000"              "tks93_50000"             
## [39] "satelitkarte_tks93_50000" "topo42_10000"            
## [41] "topo42_100000"            "topo42_200000"           
## [43] "topo42_25000"             "topo42_50000"            
## [45] "topo63_10000"             "latvija1940_75000"       
## [47] "topo63_25000"             "tks93_100000"            
## [49] "tks93_200000"             "tks93_10000"             
## [51] "ainavzemes_poly"          "ainavapvidi_poly"        
## [53] "autoceli_line"           
## attr(,"driver")
## [1] "OpenFileGDB"
## attr(,"nlayers")
## [1] 53

Mūs šajā gadījumā interesē tieši slānis “novadi_2009_pol” (novadi doti kā poligoni). Slāņa importēšanai programmā R izmanto funkciju readOGR(), norādot datubāzes atrašanās vietu un slāņa nosaukumu.

novadi <- readOGR(dsn=Latvijagdb,layer="novadi_2009_pol")

## OGR data source with driver: OpenFileGDB 
## Source: "/Users/didzis/Dropbox/DidzaDati/eksper/GIS_Latvija10.2/GIS_Latvija10.2.gdb", layer: "novadi_2009_pol"
## with 119 features
## It has 8 fields

names(novadi)

## [1] "Tips"         "Nosaukums"    "Nos_pilns"    "ATVK"         "Plan_reg"    
## [6] "skaits"       "Shape_Length" "Shape_Area"

Lai varētu veidot kartes programmā R, importētais datu slānis, kas ir kā telpiskais objekts, ir jāpārvēršs par datu tabulu (dataframe) ar funkciju fortify() no paketes ggplot2. Funkcijā fortify() ir jānorāda, kurš mainīgais ir identifikators objektiem. Šajā gadījumā tas ir “Nos_pilns”.

library(ggplot2)
novadi.df<-fortify(novadi,region="Nos_pilns")

Tagad var importēt (izveidot) datu tabulu, kurā ir tādi paši novadu nosaukumi kā kartes datos. Piemēra tabula tiek izveidota ar divām datu kolonnām - vienā skaitļi, otrā šie paši skaitļi sadalīti līmeņos.

set.seed(123)
dati <- data.frame(Novads=unique(novadi.df$id),Skaitli=runif(119,100,1000))
dati$Limeni <- cut(dati$Skaitli,breaks = 5)
head(dati)

##               Novads  Skaitli     Limeni
## 1       Ādažu novads 358.8198  (279,458]
## 2     Aglonas novads 809.4746  (637,816]
## 3 Aizkraukles novads 468.0792  (458,637]
## 4    Aizputes novads 894.7157  (816,996]
## 5    Aknīstes novads 946.4206  (816,996]
## 6      Alojas novads 141.0008 (99.7,279]

Kartes datus un informācijas tabulu nepieciešams savienot kopā ar funkciju merge() un norādot, kurās kolonnās ir novada nosaukumi katrā no tabulām.

kopa <- merge(novadi.df,dati,by.x="id",by.y="Novads")
head(kopa)

##             id     long      lat order  hole piece          group  Skaitli
## 1 Ādažu novads 530743.9 328340.0     1 FALSE     1 Ādažu novads.1 358.8198
## 2 Ādažu novads 530729.8 328333.0     2 FALSE     1 Ādažu novads.1 358.8198
## 3 Ādažu novads 530596.0 328184.4     3 FALSE     1 Ādažu novads.1 358.8198
## 4 Ādažu novads 530582.7 328114.4     4 FALSE     1 Ādažu novads.1 358.8198
## 5 Ādažu novads 530461.8 327890.0     5 FALSE     1 Ādažu novads.1 358.8198
## 6 Ādažu novads 530341.0 327820.9     6 FALSE     1 Ādažu novads.1 358.8198
##      Limeni
## 1 (279,458]
## 2 (279,458]
## 3 (279,458]
## 4 (279,458]
## 5 (279,458]
## 6 (279,458]

Datu attēlošanai izmanto geom_polygon(), kurā aizpildījums ir atkarīgs no datu kolonnas. Atkarībā no tā, vai datu kolonna ir skaitliska, vai kategorijas, mainās izskats un leģenda. Attiecībā uz Latvijas novadu karti, geom_polygon() ir jāizmanto divas reizes, jo Daugavpils un Rēzeknes pilsētas ir aptvertas ar atbilstošajiem novadiem, tāpēc, lai šo pilsētu dati nepazustu, tie jāūzliek “pa virsu” atsevišķi. Tā kā šajos datos ir izmantot LKS koordinātes, tad varam izmanto coord_fixed(), lai panāktu, ka karte tiek attēlota proporcionāli. Ar theme_void() uzreiz var noņemt visas asis un citu nevajadzīgo noformējumu.

ggplot(kopa,aes(long,lat,group=group))+
      geom_polygon(aes(fill=Skaitli),color="white") + 
      geom_polygon(data=subset(kopa,id %in% c("Rēzekne","Daugavpils")),
                   aes(fill=Skaitli),color="white") +
      scale_fill_gradient(low = "lightgreen", high = "darkgreen" )+
      coord_fixed() + theme_void()

ggplot(kopa,aes(long,lat,group=group))+
      geom_polygon(aes(fill=Limeni),color="black",size=0.2) + 
      geom_polygon(data=subset(kopa,id %in% c("Rēzekne","Daugavpils")),aes(fill=Limeni),color="black",size=0.2) + 
      scale_fill_brewer()+
      coord_fixed() + theme_void()