library(maptools)
## Loading required package: sp
## Checking rgeos availability: TRUE
library(sp)
library(RColorBrewer)
library(classInt)
Lectura del shapefile
muni.map <- readShapePoly(fn = "~/Descargas/DERA_actualizado/G17_Division_Administrativa/da02_term_munic.shp")
# vemos que tiene
slotNames(muni.map)
## [1] "data" "polygons" "plotOrder" "bbox" "proj4string"
head(slot(muni.map, "data"))
## COD_MUN MUNICIPIO PROVINCIA COD_ENT
## 0 04001 Abla Almer\xeda d06
## 1 04002 Abrucena Almer\xeda d06
## 2 04003 Adra Almer\xeda d06
## 3 04004 Alb\xe1nchez Almer\xeda d06
## 4 04005 Alboloduy Almer\xeda d06
## 5 04006 Albox Almer\xeda d06
# Ponemos encoding a latin1
Encoding(levels(muni.map$MUNICIPIO)) <- "latin1"
Encoding(levels(muni.map$PROVINCIA)) <- "latin1"
head(slot(muni.map, "data"))
## COD_MUN MUNICIPIO PROVINCIA COD_ENT
## 0 04001 Abla Almería d06
## 1 04002 Abrucena Almería d06
## 2 04003 Adra Almería d06
## 3 04004 Albánchez Almería d06
## 4 04005 Alboloduy Almería d06
## 5 04006 Albox Almería d06
Unimos los polígonos con los códigos de los municipios (si los uno con los de las provincias, obtengo mapa provincial) de esta forma ahora el ID de los polígonos es el código del municipio
anda.map <- unionSpatialPolygons(muni.map, muni.map$COD_MUN)
## Loading required package: rgeos
## rgeos version: 0.3-3, (SVN revision 437)
## GEOS runtime version: 3.3.3-CAPI-1.7.4
## Polygon checking: TRUE
Para el mapa sin unir el ID de los polígonos va de 0 a número de polígonos
muni.map@polygons[[1]]@ID
## [1] "0"
muni.map@polygons[[2]]@ID
## [1] "1"
Al unir con unionSpatialPolygons les ponemos como ID el COD_MUN
anda.map@polygons[[1]]@ID
## [1] "04001"
Puedo extraer el ID para todos los polígonos
# con sapply aplico una función a cada uno de los polígonos. la función que
# aplico es la de extraer el ID sppply, simplifica la salida si es posible.
# es decir si puede devolver un vector en vez de una lista, lo hará
sapply(anda.map@polygons, function(x) x@ID)
## [1] "04001" "04002" "04003" "04004" "04005" "04006" "04007" "04008"
## [9] "04009" "04010" "04011" "04012" "04013" "04014" "04015" "04016"
## [17] "04017" "04018" "04019" "04020" "04021" "04022" "04023" "04024"
## [25] "04026" "04027" "04028" "04029" "04030" "04031" "04032" "04033"
## [33] "04034" "04035" "04036" "04037" "04038" "04041" "04043" "04044"
## [41] "04045" "04046" "04047" "04048" "04049" "04050" "04051" "04052"
## [49] "04053" "04054" "04055" "04056" "04057" "04058" "04059" "04060"
## [57] "04061" "04062" "04063" "04064" "04065" "04066" "04067" "04068"
## [65] "04069" "04070" "04071" "04072" "04073" "04074" "04075" "04076"
## [73] "04077" "04078" "04079" "04080" "04081" "04082" "04083" "04084"
## [81] "04085" "04086" "04087" "04088" "04089" "04090" "04091" "04092"
## [89] "04093" "04094" "04095" "04096" "04097" "04098" "04099" "04100"
## [97] "04101" "04102" "04103" "04901" "04902" "04903" "11001" "11002"
## [105] "11003" "11004" "11005" "11006" "11007" "11008" "11009" "11010"
## [113] "11011" "11012" "11013" "11014" "11015" "11016" "11017" "11018"
## [121] "11019" "11020" "11021" "11022" "11023" "11024" "11025" "11026"
## [129] "11027" "11028" "11029" "11030" "11031" "11032" "11033" "11034"
## [137] "11035" "11036" "11037" "11038" "11039" "11040" "11041" "11042"
## [145] "11901" "11902" "14001" "14002" "14003" "14004" "14005" "14006"
## [153] "14007" "14008" "14009" "14010" "14011" "14012" "14013" "14014"
## [161] "14015" "14016" "14017" "14018" "14019" "14020" "14021" "14022"
## [169] "14023" "14024" "14025" "14026" "14027" "14028" "14029" "14030"
## [177] "14031" "14032" "14033" "14034" "14035" "14036" "14037" "14038"
## [185] "14039" "14040" "14041" "14042" "14043" "14044" "14045" "14046"
## [193] "14047" "14048" "14049" "14050" "14051" "14052" "14053" "14054"
## [201] "14055" "14056" "14057" "14058" "14059" "14060" "14061" "14062"
## [209] "14063" "14064" "14065" "14066" "14067" "14068" "14069" "14070"
## [217] "14071" "14072" "14073" "14074" "14075" "18001" "18002" "18003"
## [225] "18004" "18005" "18006" "18007" "18010" "18011" "18012" "18013"
## [233] "18014" "18015" "18016" "18017" "18018" "18020" "18021" "18022"
## [241] "18023" "18024" "18025" "18027" "18028" "18029" "18030" "18032"
## [249] "18033" "18034" "18035" "18036" "18037" "18038" "18039" "18040"
## [257] "18042" "18043" "18044" "18045" "18046" "18047" "18048" "18049"
## [265] "18050" "18051" "18053" "18054" "18056" "18057" "18059" "18061"
## [273] "18062" "18063" "18064" "18066" "18067" "18068" "18069" "18070"
## [281] "18071" "18072" "18074" "18076" "18078" "18079" "18082" "18083"
## [289] "18084" "18085" "18086" "18087" "18088" "18089" "18093" "18094"
## [297] "18095" "18096" "18097" "18098" "18099" "18100" "18101" "18102"
## [305] "18103" "18105" "18107" "18108" "18109" "18111" "18112" "18114"
## [313] "18115" "18116" "18117" "18119" "18120" "18121" "18122" "18123"
## [321] "18124" "18126" "18127" "18128" "18132" "18133" "18134" "18135"
## [329] "18136" "18137" "18138" "18140" "18141" "18143" "18144" "18145"
## [337] "18146" "18147" "18148" "18149" "18150" "18151" "18152" "18153"
## [345] "18154" "18157" "18158" "18159" "18161" "18162" "18163" "18164"
## [353] "18165" "18167" "18168" "18170" "18171" "18173" "18174" "18175"
## [361] "18176" "18177" "18178" "18179" "18180" "18181" "18182" "18183"
## [369] "18184" "18185" "18187" "18188" "18189" "18192" "18193" "18194"
## [377] "18901" "18902" "18903" "18904" "18905" "18906" "18907" "18908"
## [385] "18909" "18910" "18911" "18912" "18913" "21001" "21002" "21003"
## [393] "21004" "21005" "21006" "21007" "21008" "21009" "21010" "21011"
## [401] "21012" "21013" "21014" "21015" "21016" "21017" "21018" "21019"
## [409] "21020" "21021" "21022" "21023" "21024" "21025" "21026" "21027"
## [417] "21028" "21029" "21030" "21031" "21032" "21033" "21034" "21035"
## [425] "21036" "21037" "21038" "21039" "21040" "21041" "21042" "21043"
## [433] "21044" "21045" "21046" "21047" "21048" "21049" "21050" "21051"
## [441] "21052" "21053" "21054" "21055" "21056" "21057" "21058" "21059"
## [449] "21060" "21061" "21062" "21063" "21064" "21065" "21066" "21067"
## [457] "21068" "21069" "21070" "21071" "21072" "21073" "21074" "21075"
## [465] "21076" "21077" "21078" "21079" "23001" "23002" "23003" "23004"
## [473] "23005" "23006" "23007" "23008" "23009" "23010" "23011" "23012"
## [481] "23014" "23015" "23016" "23017" "23018" "23019" "23020" "23021"
## [489] "23024" "23025" "23026" "23027" "23028" "23029" "23030" "23031"
## [497] "23032" "23033" "23034" "23035" "23037" "23038" "23039" "23040"
## [505] "23041" "23042" "23043" "23044" "23045" "23046" "23047" "23048"
## [513] "23049" "23050" "23051" "23052" "23053" "23054" "23055" "23056"
## [521] "23057" "23058" "23059" "23060" "23061" "23062" "23063" "23064"
## [529] "23065" "23066" "23067" "23069" "23070" "23071" "23072" "23073"
## [537] "23074" "23075" "23076" "23077" "23079" "23080" "23081" "23082"
## [545] "23084" "23085" "23086" "23087" "23088" "23090" "23091" "23092"
## [553] "23093" "23094" "23095" "23096" "23097" "23098" "23099" "23101"
## [561] "23901" "23902" "23903" "23904" "23905" "29001" "29002" "29003"
## [569] "29004" "29005" "29006" "29007" "29008" "29009" "29010" "29011"
## [577] "29012" "29013" "29014" "29015" "29016" "29017" "29018" "29019"
## [585] "29020" "29021" "29022" "29023" "29024" "29025" "29026" "29027"
## [593] "29028" "29029" "29030" "29031" "29032" "29033" "29034" "29035"
## [601] "29036" "29037" "29038" "29039" "29040" "29041" "29042" "29043"
## [609] "29044" "29045" "29046" "29047" "29048" "29049" "29050" "29051"
## [617] "29052" "29053" "29054" "29055" "29056" "29057" "29058" "29059"
## [625] "29060" "29061" "29062" "29063" "29064" "29065" "29066" "29067"
## [633] "29068" "29069" "29070" "29071" "29072" "29073" "29074" "29075"
## [641] "29076" "29077" "29079" "29080" "29081" "29082" "29083" "29084"
## [649] "29085" "29086" "29087" "29088" "29089" "29090" "29091" "29092"
## [657] "29093" "29094" "29095" "29096" "29097" "29098" "29099" "29100"
## [665] "29901" "29902" "41001" "41002" "41003" "41004" "41005" "41006"
## [673] "41007" "41008" "41009" "41010" "41011" "41012" "41013" "41014"
## [681] "41015" "41016" "41017" "41018" "41019" "41020" "41021" "41022"
## [689] "41023" "41024" "41025" "41026" "41027" "41028" "41029" "41030"
## [697] "41031" "41032" "41033" "41034" "41035" "41036" "41037" "41038"
## [705] "41039" "41040" "41041" "41042" "41043" "41044" "41045" "41046"
## [713] "41047" "41048" "41049" "41050" "41051" "41052" "41053" "41054"
## [721] "41055" "41056" "41057" "41058" "41059" "41060" "41061" "41062"
## [729] "41063" "41064" "41065" "41066" "41067" "41068" "41069" "41070"
## [737] "41071" "41072" "41073" "41074" "41075" "41076" "41077" "41078"
## [745] "41079" "41080" "41081" "41082" "41083" "41084" "41085" "41086"
## [753] "41087" "41088" "41089" "41090" "41091" "41092" "41093" "41094"
## [761] "41095" "41096" "41097" "41098" "41099" "41100" "41101" "41102"
## [769] "41901" "41902" "41903"
Ahora puedo seleccionar polígonos por su código y dibujarlos Dibujo el término municipal de mi pueblo, que se corresponde con el código 18053. O dibujar polígonos sueltos
plot(anda.map["18053"]) # mi pueblo
plot(anda.map[c("18053", "18050")]) # mi pueblo y otro
O ponerlos de un color
plot(anda.map)
plot(anda.map["18053"], add = TRUE, col = "lightblue")
¿Como coloreo si tengo unos datos externos al mapa? A partir de un csv que tenga los datos
datos <- read.csv("~/Dropbox/Maps_and_mrp_article/Mapas_en_R/mapas_OPAM/padron/datos/municipios/dat_map_1.csv")
head(datos)
## muni cod_mun tot.pob tot.ext
## 1 04001-Abla 4001 1465 177
## 2 04002-Abrucena 4002 1360 109
## 3 04003-Adra 4003 24626 3063
## 4 04004-Albánchez 4004 774 281
## 5 04005-Alboloduy 4005 674 15
## 6 04006-Albox 4006 10821 2741
Si los datos están en el mismo orden que los polígonos es fácil dibujar el mapa
Voy a dibujar el porcentaje de extranjeros Primero calculo el porcentaje
datos$porcentaje <- with(datos, 100 * tot.ext/tot.pob)
Asigno colores y divido porcentaje en 5 intervalos
# los colores lso pongo en hexadecimal, pero se podría elegir 'red','blue',
# etc
azul_1 <- "#E1E8FF"
azul_2 <- "#BAC9FE"
azul_3 <- "#7284FE"
azul_4 <- "#0228DF"
azul_5 <- "#01167C"
colores <- c(azul_1, azul_2, azul_3, azul_4, azul_5)
Divido la variable porcentaje, usando cut por ejemplo. También se pueden usar funciones de la librería classInt
corte.factor <- cut(datos$porcentaje, 5)
table(corte.factor)
## corte.factor
## (-0.0683,13.6] (13.6,27.3] (27.3,41] (41,54.7] (54.7,68.4]
## 647 82 26 11 5
Me creo un vector de longitud, el número de datos pero que contenga los valores de los colores según los cortes
colores2 <- colores[corte.factor]
head(colores2)
## [1] "#E1E8FF" "#E1E8FF" "#E1E8FF" "#7284FE" "#E1E8FF" "#BAC9FE"
Como mis datos están en el mismo orden que los polígonos dibujar el mapa es tán fácil como
plot(anda.map, col = colores2)
Si no estuvieran en el mismo orden, podría usar merge para unir los colores a muni.map@data (muni.map es el mapa antes de hacer el unionSpatialPolygons) utilizando una variable común, como el código del municipio. ¡¡OJO!! que al cargar los datos el COD_MUN lo ha considerado numérico en vez de como cáracter,y los códigos que empezaban por 0 como 04006, los toma como 4006.
Podemos solucionarlo utilizando colClasses al leer el csv
datos <- read.csv("~/Dropbox/Maps_and_mrp_article/Mapas_en_R/mapas_OPAM/padron/datos/municipios/dat_map_1.csv",
colClasses = c(rep("factor", 2), rep("integer", 2)))
datos$porcentaje <- with(datos, 100 * tot.ext/tot.pob)
corte.factor <- cut(datos$porcentaje, 5)
datos$colores2 <- colores[corte.factor]
head(datos$colores2)
## [1] "#E1E8FF" "#E1E8FF" "#E1E8FF" "#7284FE" "#E1E8FF" "#BAC9FE"
Veo el nombre de muni.map@data (esto equivale al dbf de la capa)
names(muni.map@data)
## [1] "COD_MUN" "MUNICIPIO" "PROVINCIA" "COD_ENT"
le uno el data.frame datos
muni.map@data <- merge(muni.map@data, datos, by.x = "COD_MUN", by.y = "cod_mun")
head(muni.map@data)
## COD_MUN MUNICIPIO PROVINCIA COD_ENT muni tot.pob tot.ext
## 1 04001 Abla Almería d06 04001-Abla 1465 177
## 2 04002 Abrucena Almería d06 04002-Abrucena 1360 109
## 3 04003 Adra Almería d06 04003-Adra 24626 3063
## 4 04004 Albánchez Almería d06 04004-Albánchez 774 281
## 5 04005 Alboloduy Almería d06 04005-Alboloduy 674 15
## 6 04006 Albox Almería d06 04006-Albox 10821 2741
## porcentaje colores2
## 1 12.082 #E1E8FF
## 2 8.015 #E1E8FF
## 3 12.438 #E1E8FF
## 4 36.305 #7284FE
## 5 2.226 #E1E8FF
## 6 25.330 #BAC9FE
names(muni.map) # ya está colores2
## [1] "COD_MUN" "MUNICIPIO" "PROVINCIA" "COD_ENT" "muni"
## [6] "tot.pob" "tot.ext" "porcentaje" "colores2"
Dibujamos
plot(muni.map, col = muni.map$colores2)