This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
En el siguiente trabajo se realizará la lectura de un shepfile descargado de la base de datos del DANE del geoportal, de Colombia y un filtrado del departamento de Cundinamarca
#install.packages("tidyverse")
#install.packages("sf")
library(tidyverse)
library(sf)
Para realizar esta lectura de datos lo primero que se debe hacer es descargar la base de datos con la que se va a trabajar, para este caso los datos fueron descargados de la base de datos del geoportal del DANE, se crea una variable llamada “deptos” y utilizamos la función read_sf("") de la librería previamente instalada (sf) y añadimos la ubicación del archivo con la base de datos
deptos <- read_sf("C:/Users/mario/Documents/Ing. Agronomica/9 Semestre/Geomatica/Cundinamarca/COL_adm/COL_adm1.shp")
Ahora utilizando la función head() llamamos la tabla de datos del archivo descargado
head(deptos)
## Simple feature collection with 6 features and 9 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -77.149 ymin: -4.228429 xmax: -69.36835 ymax: 11.10792
## CRS: 4326
## # A tibble: 6 x 10
## ID_0 ISO NAME_0 ID_1 NAME_1 TYPE_1 ENGTYPE_1 NL_NAME_1 VARNAME_1
## <dbl> <chr> <chr> <dbl> <chr> <chr> <chr> <chr> <chr>
## 1 53 COL Colom~ 1 Amazo~ Comis~ Commissi~ <NA> <NA>
## 2 53 COL Colom~ 2 Antio~ Depar~ Departme~ <NA> <NA>
## 3 53 COL Colom~ 3 Arauca Inten~ Intendan~ <NA> <NA>
## 4 53 COL Colom~ 4 Atlán~ Depar~ Departme~ <NA> <NA>
## 5 53 COL Colom~ 5 Bolív~ Depar~ Departme~ <NA> <NA>
## 6 53 COL Colom~ 6 Boyacá Depar~ Departme~ <NA> <NA>
## # ... with 1 more variable: geometry <MULTIPOLYGON [°]>
Para consultar el sistema de referencia de coordenadas, se utiliza la siguiente función:
st_crs(deptos)
## Coordinate Reference System:
## User input: 4326
## wkt:
## GEOGCS["GCS_WGS_1984",
## DATUM["WGS_1984",
## SPHEROID["WGS_84",6378137.0,298.257223563]],
## PRIMEM["Greenwich",0.0],
## UNIT["Degree",0.0174532925199433],
## AUTHORITY["EPSG","4326"]]
Este sistema de referencia de coordenadas, utiliza el datum del Sistema geodésico mundial 1984.
Para realizar esta visualización utilizamos el siguiente código:
ggplot() + geom_sf(data = deptos)
Para realizar una proyección con el sistema de coordenadas EPSG:3978, utilizamos el siguiente código
ggplot() + geom_sf(data = deptos) + coord_sf(crs=st_crs(3978))
El anterior mapa tiene un sistema de coordenadas proyectado EPSG:3978 NAD83 / Canada Atlas Lambert, por esto realizaremos una transformación del sistema de coordenadas EPSG:3978 a EPSG:32618WGS 84 / UTM zone 18N, utilizando el siguiente código:
deptos_utm <- st_transform(deptos, crs = st_crs(32618))
Ahora llamamos a los datos (deptos_utm) que ya estan transformados al sistema de coordenadas EPSG:32618 / UTM 18N
deptos_utm
## Simple feature collection with 32 features and 9 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -245935.3 ymin: -469204.3 xmax: 1407491 ymax: 1763314
## CRS: EPSG:32618
## # A tibble: 32 x 10
## ID_0 ISO NAME_0 ID_1 NAME_1 TYPE_1 ENGTYPE_1 NL_NAME_1 VARNAME_1
## <dbl> <chr> <chr> <dbl> <chr> <chr> <chr> <chr> <chr>
## 1 53 COL Colom~ 1 Amazo~ Comis~ Commissi~ <NA> <NA>
## 2 53 COL Colom~ 2 Antio~ Depar~ Departme~ <NA> <NA>
## 3 53 COL Colom~ 3 Arauca Inten~ Intendan~ <NA> <NA>
## 4 53 COL Colom~ 4 Atlán~ Depar~ Departme~ <NA> <NA>
## 5 53 COL Colom~ 5 Bolív~ Depar~ Departme~ <NA> <NA>
## 6 53 COL Colom~ 6 Boyacá Depar~ Departme~ <NA> <NA>
## 7 53 COL Colom~ 7 Córdo~ Depar~ Departme~ <NA> <NA>
## 8 53 COL Colom~ 8 Caldas Depar~ Departme~ <NA> <NA>
## 9 53 COL Colom~ 9 Caque~ Inten~ Intendan~ <NA> <NA>
## 10 53 COL Colom~ 10 Casan~ Inten~ Intendan~ <NA> <NA>
## # ... with 22 more rows, and 1 more variable: geometry <MULTIPOLYGON [m]>
A continuación realizaremos el mapa usando el código anteriormente usado pero en este caso se utilizaran los datos transformados
ggplot() + geom_sf(data = deptos_utm)
Como solo estamos interesados en los datos del departamento de Cundinamarca podemos filtrar la base de datos: Para esto crearemos una variable llamada “Cundinamarca” Y utilizaremos el siguiente código “deptos %>% filter(NAME_1 ==”Cundinamarca“), donde el filtro se hace en la base de datos donde estan los departamentos, como se observa en los datos anteriores el encabezado de los departamentos es”NAME_1"
cundinamarca <- deptos %>% filter(NAME_1 == "Cundinamarca")
Ahora realizamos el mapa del departamento filtrado, utilizando el siguiente código
ggplot() + geom_sf(data = cundinamarca)
Añadimos la base de dats de todos los departamentos
munic <- read_sf("C:/Users/mario/Documents/Ing. Agronomica/9 Semestre/Geomatica/Cundinamarca/COL_adm/COL_adm2.shp")
munic
## Simple feature collection with 1065 features and 11 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -81.84153 ymin: -4.228429 xmax: -66.87033 ymax: 15.91247
## CRS: 4326
## # A tibble: 1,065 x 12
## ID_0 ISO NAME_0 ID_1 NAME_1 ID_2 NAME_2 TYPE_2 ENGTYPE_2 NL_NAME_2
## <dbl> <chr> <chr> <dbl> <chr> <dbl> <chr> <chr> <chr> <chr>
## 1 53 COL Colom~ 1 Amazo~ 1 El En~ Corre~ Corregim~ <NA>
## 2 53 COL Colom~ 1 Amazo~ 2 La Ch~ Corre~ Corregim~ <NA>
## 3 53 COL Colom~ 1 Amazo~ 3 La Pe~ Corre~ Corregim~ <NA>
## 4 53 COL Colom~ 1 Amazo~ 4 Letic~ Munic~ Municipa~ <NA>
## 5 53 COL Colom~ 1 Amazo~ 5 Mirit~ Corre~ Corregim~ <NA>
## 6 53 COL Colom~ 1 Amazo~ 6 Puert~ Munic~ Municipa~ <NA>
## 7 53 COL Colom~ 1 Amazo~ 7 Puert~ Corre~ Corregim~ <NA>
## 8 53 COL Colom~ 1 Amazo~ 8 Tarap~ Corre~ Corregim~ <NA>
## 9 53 COL Colom~ 2 Antio~ 9 Abejo~ Munic~ Municipa~ <NA>
## 10 53 COL Colom~ 2 Antio~ 10 Abria~ Munic~ Municipa~ <NA>
## # ... with 1,055 more rows, and 2 more variables: VARNAME_2 <chr>,
## # geometry <MULTIPOLYGON [°]>
Filtramos los datos de solamente el departamento de Cundinamarca
mun_cundinamarca <- munic %>% filter(NAME_1 == "Cundinamarca")
ggplot() + geom_sf(data = mun_cundinamarca)
mun_cundinamarca
## Simple feature collection with 115 features and 11 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -74.919 ymin: 3.6546 xmax: -73.0677 ymax: 5.820602
## CRS: 4326
## # A tibble: 115 x 12
## ID_0 ISO NAME_0 ID_1 NAME_1 ID_2 NAME_2 TYPE_2 ENGTYPE_2 NL_NAME_2
## * <dbl> <chr> <chr> <dbl> <chr> <dbl> <chr> <chr> <chr> <chr>
## 1 53 COL Colom~ 14 Cundi~ 490 Agua ~ Munic~ Municipa~ <NA>
## 2 53 COL Colom~ 14 Cundi~ 491 Albán Munic~ Municipa~ <NA>
## 3 53 COL Colom~ 14 Cundi~ 492 Anapo~ Munic~ Municipa~ <NA>
## 4 53 COL Colom~ 14 Cundi~ 493 Anola~ Munic~ Municipa~ <NA>
## 5 53 COL Colom~ 14 Cundi~ 494 Apulo Munic~ Municipa~ <NA>
## 6 53 COL Colom~ 14 Cundi~ 495 Arbel~ Munic~ Municipa~ <NA>
## 7 53 COL Colom~ 14 Cundi~ 496 Beltr~ Munic~ Municipa~ <NA>
## 8 53 COL Colom~ 14 Cundi~ 497 Bitui~ Munic~ Municipa~ <NA>
## 9 53 COL Colom~ 14 Cundi~ 498 Bojacá Munic~ Municipa~ <NA>
## 10 53 COL Colom~ 14 Cundi~ 499 Cáque~ Munic~ Municipa~ <NA>
## # ... with 105 more rows, and 2 more variables: VARNAME_2 <chr>,
## # geometry <MULTIPOLYGON [°]>
cundinamarca_points<- st_centroid(mun_cundinamarca)
cundinamarca_points <- cbind(mun_cundinamarca, st_coordinates(st_centroid(mun_cundinamarca$geometry)))
ggplot(cundinamarca) +
geom_sf() +
geom_sf(data = cundinamarca_points, fill = "antiquewhite") +
geom_text(data = cundinamarca_points, aes(x=X, y=Y,label = ID_2), size = 2, ) +
coord_sf(xlim = c(-75.5, -72.5), ylim = c(3.5, 6), expand = FALSE)
library(scales)
ggplot(cundinamarca) +
geom_sf(data=cundinamarca_points, aes(x=X, y=Y, fill =
ID_2), color = "black", size = 0.25) +
geom_text(data = cundinamarca_points, aes(x=X, y=Y,label = ID_2), size = 2) +
theme(aspect.ratio=1)+ scale_fill_continuous(name="ID_2",breaks = pretty_breaks(n = 5))+ labs(title="Cundinamarca")
ggsave("cundinamarca_municipios.pdf")
ggsave("map_cundinamarca.png", width = 6, height = 6, dpi = "screen")
#install.packages("leaflet")
library(leaflet)
ant_points <- as(cundinamarca_points, 'Spatial')
#install.packages("lwgeom")
library(lwgeom)
mun_cundinamarca$area <- st_area(mun_cundinamarca) #Take care of units
mun_cundinamarca$km2 <- mun_cundinamarca$area/(1000000)
ant_mun <- as(mun_cundinamarca, 'Spatial')
##ant_mun
bins <- c(0, 50, 100, 200, 300, 500, 1000, 2000, Inf)
labels <- mun_cundinamarca$NAME_2
labels
## [1] "Agua de Dios" "Albán"
## [3] "Anapoima" "Anolaima"
## [5] "Apulo" "Arbeláez"
## [7] "Beltrán" "Bituima"
## [9] "Bojacá" "Cáqueza"
## [11] "Cabrera" "Cachipay"
## [13] "Cajicá" "Caparrapí"
## [15] "Carmen de Carupa" "Chía"
## [17] "Chaguaní" "Chipaque"
## [19] "Choachí" "Chocontá"
## [21] "Cogua" "Cota"
## [23] "Cucunubá" "El Colegio"
## [25] "El Peñon" "Fómeque"
## [27] "Facatativá" "Fúquene"
## [29] "Fosca" "Funza"
## [31] "Fusagasugá" "Gachalá"
## [33] "Gachancipá" "Gachetá"
## [35] "Gama" "Girardot"
## [37] "Guachetá" "Guaduas"
## [39] "Guasca" "Guataquí"
## [41] "Guatavita" "Guayabal de Síquima"
## [43] "Guayabetal" "Gutiérrez"
## [45] "Jerusalén" "Junín"
## [47] "La Calera" "La Mesa"
## [49] "La Palma" "La Peña"
## [51] "La Vega" "Lenguazaque"
## [53] "Machetá" "Madrid"
## [55] "Manta" "Medina"
## [57] "Mosquera" "Nariño"
## [59] "Nemocón" "Nilo"
## [61] "Nimaima" "Nocaima"
## [63] "Pacho" "Paime"
## [65] "Pandi" "Paratebueno"
## [67] "Pasca" "Puerto Salgar"
## [69] "Pulí" "Quebradanegra"
## [71] "Quetame" "Quipile"
## [73] "Ricaurte" "San Antonio del Tequendama"
## [75] "San Bernardo" "San Cayetano"
## [77] "San Francisco" "San Juan de Río Seco"
## [79] "Santafé de Bogotá" "Sasaima"
## [81] "Sesquilé" "Sibaté"
## [83] "Silvania" "Simijaca"
## [85] "Soacha" "Sopó"
## [87] "Subachoque" "Suesca"
## [89] "Supatá" "Susa"
## [91] "Sutatausa" "Tabio"
## [93] "Tausa" "Tena"
## [95] "Tenjo" "Tibacuy"
## [97] "Tibirita" "Tocaima"
## [99] "Tocancipá" "Topaipí"
## [101] "Ubalá" "Ubaque"
## [103] "Une" "Utica"
## [105] "Venecia" "Vergara"
## [107] "Vianí" "Villa de San Diego de Ubaté"
## [109] "Villagómez" "Villapinzón"
## [111] "Villeta" "Viotá"
## [113] "Yacopí" "Zipacón"
## [115] "Zipaquirá"
library(leaflet)
library(lwgeom)
b <- leaflet(ant_mun) %>%
setView(-74, 4.7, 7.5) %>% addPolygons (data = ant_mun,
weight = 2,
opacity = 1,
color = "green",
dashArray = "3",
fillOpacity = 0.7,
highlight = highlightOptions(
weight = 5,
color = "#666",
dashArray = "",
fillOpacity = 0.7,
bringToFront = TRUE))
b
leaflet() %>%
addProviderTiles(providers$Esri.WorldImagery, options= providerTileOptions(opacity = 0.99)) %>%
addPolygons(data = ant_mun, popup= ant_mun$NAME_2,
stroke = TRUE, fillOpacity = 0.25, smoothFactor = 0.25
)