MAPA TEMÁTICO DE TOLIMA
Valentina Diaz Diaz
10/10/2023
1. Introducción
Este es el segundo cuaderno R Markdown para el curso Geomatica Básica 2022. Ilustra cómo hacer mapas temáticos que muestren la participación municipal de los dos grupos de cultivos más importantes para un departamento determinado. Usaremos como fuente principal los archivos csv guardados en el cuaderno EVA, así como un shapefile de municipios obtenidos en clase.
2. Configuración
En este paso, instalaremos y cargaremos las bibliotecas R necesarias.
#run the following lines from the command line:
#install.packages('dplyr')
#install.packages('readxl')
#install.packages('sf)library(sf)## Linking to GEOS 3.11.2, GDAL 3.6.2, PROJ 9.2.0; sf_use_s2() is TRUElibrary(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(readr)3. Leer archivos relacionados con municipios, cultivos y ciudades
Trabajaré con el Departamento de Tolima.
Antes de comenzar, asegúrese de haber guardado los siguientes archivos en su computadora:
-Un archivo de forma con los municipios de su departamento (como se produjo en clase) -Un archivo csv con estadísticas de EVA 2020 para el grupo de cultivos seleccionado (como se produjo en el primer cuaderno) Además, se necesitara un archivo con ciudades de Colombia. Puede descargarlo desde https://simplemaps.com/data/co-cities. Luego, mueva el archivo co.csv descargado a su directorio de trabajo.
-Ahora, procedamos a leer los archivos basados en EVA obtenidos del primer cuaderno:
library(readxl)
cereales <- read_excel("EVA_19_22.xlsx")
View(cereales)alg_tol=cereales%>%group_by(COD_MUN, MUNICIPIO, GRUPO)%>%filter(DEPTO=="Tolima"& GRUPO=="Cereales")%>%
summarize(max_prod = max(PRODUCCION, na.rm = TRUE)) %>%
arrange(desc(max_prod))## `summarise()` has grouped output by 'COD_MUN', 'MUNICIPIO'. You can override
## using the `.groups` argument.write_csv(alg_tol,"Tolima_Algodon.csv")- Lee los municipios de tu departamento:
library(sf)# this shapefile must have a 4326 EPSG code
# this shapefile was obtained in class using QGIS
(mun.tmp = st_read("C:/Users/camil/OneDrive/Escritorio/Geomatica/CuadernoEVA/MGN_MPIO_POLITICO.shp"))## Reading layer `MGN_MPIO_POLITICO' from data source
## `C:\Users\camil\OneDrive\Escritorio\Geomatica\CuadernoEVA\MGN_MPIO_POLITICO.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 1121 features and 12 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -81.73562 ymin: -4.229406 xmax: -66.84722 ymax: 13.39473
## Geodetic CRS: MAGNA-SIRGAS-Seleccionemos algunos atributos para limpiar el objeto:
# check the objet output in the last chunk and
# change attribute names according to your own data
mun.tmp %>% filter(DPTO_CNMBR=="TOLIMA") %>% select(MPIO_CCDGO, MPIO_CNMBR, MPIO_NAREA) -> municipiosmunicipiosst_transform(x = municipios,crs = 4326)-Ahora, lea el archivo csv de las ciudades:
ciudadesco=read.csv("ciudadesco.csv")ciudadescoTenga en cuenta que hay 1102 filas que representan ciudades en Colombia. Tenga en cuenta que el objeto de ciudades es un objeto no espacial (un tibble, en realidad).
Tenga en cuenta también que los valores de las coordenadas se almacenan en los atributos “lng” (latitud) y “lat” (longitud).
Para convertir ciudades en un objeto espacial (un objeto de característica simple, en este caso):
# Convert data frame to simple feature object
sf.cities <- st_as_sf(x = ciudadesco,
coords = c("lng", "lat"))sf.cities¡Tenga en cuenta el CRS que falta!
Agreguemos el CRS usando el código EPSG:
st_crs(sf.cities) <- 4326
st_crs(municipios)=4326## Warning: st_crs<- : replacing crs does not reproject data; use st_transform for
## that4. Subconjunto de datos relevantes para nuestro departamento
Como estamos interesados en un solo departamento, necesitamos crear una unión espacial:
# find points (cities) within polygons (our municipalities)
sf.cities.joined <- st_join(sf.cities, municipios, join = st_intersects)
tol.cities=filter(sf.cities.joined,admin_name=="Tolima")library(tmap)## The legacy packages maptools, rgdal, and rgeos, underpinning the sp package,
## which was just loaded, were retired in October 2023.
## Please refer to R-spatial evolution reports for details, especially
## https://r-spatial.org/r/2023/05/15/evolution4.html.
## It may be desirable to make the sf package available;
## package maintainers should consider adding sf to Suggests:.## Breaking News: tmap 3.x is retiring. Please test v4, e.g. with
## remotes::install_github('r-tmap/tmap')library(classInt)
library(ggrepel)## Loading required package: ggplot2class(alg_tol$COD_MUN)## [1] "character"class(mun.tmp$MPIO_CDPMP)## [1] "character"munic_cereales=left_join(mun.tmp,alg_tol,by=c("MPIO_CDPMP"="COD_MUN"))
munic_cereales=munic_cereales %>% filter(DPTO_CNMBR=="TOLIMA") %>% na.omit()breaks <- classIntervals(munic_cereales$max_prod, n = 6, style = 'fisher')#label breaks
lab_vec <- vector(length = length(breaks$brks)-1)
rounded_breaks <- round(breaks$brks,2)
lab_vec[1] <- paste0('[', rounded_breaks[1],' - ', rounded_breaks[2],']')
for(i in 2:(length(breaks$brks) - 1)){
lab_vec[i] <- paste0('(',rounded_breaks[i], ' - ', rounded_breaks[i+1], ']')}munic_cereales <- munic_cereales %>%
mutate(faktor_class = factor(cut(max_prod, breaks$brks, include.lowest = T), labels = lab_vec))
# change attribute name
munic_cereales$Produccion = munic_cereales$faktor_class# This code creates a new field ("mid") with centroid coordinates
munic_cereales$mid <- sf::st_centroid(munic_cereales$geometry)
# Get the longitude values
LONG = st_coordinates(munic_cereales$mid)[,1]
# Get the latitude values
LAT = st_coordinates(munic_cereales$mid)[,2]5. Haz un mapa para el grupo de cultivos más importante
ggplot(data = munic_cereales) +
geom_sf(aes(fill = Produccion)) +
geom_label_repel(aes(x = LONG, y = LAT, label = MPIO_CNMBR),
label.padding = unit(0.05,"lines"),
label.r = unit(0.025, "lines"),
label.size = 0.05) +
geom_sf(aes(fill = Produccion)) +
geom_label_repel(aes(x = LONG, y = LAT, label = MPIO_CNMBR),
label.padding = unit(0.05,"lines"),
label.r = unit(0.025, "lines"),
label.size = 0.05)
6. Haz un nuevo mapa para el segundo grupo más grande de cultivos.
-Otra opción es producir un mapa interactivo:
# see example at https://r-tmap.github.io/tmap-book/layout.html
facet = "max_prod"
cereal_map =
tm_shape(munic_cereales) + tm_polygons(facet) + tm_text(text = "MPIO_CNMBR", size = 0.7, fontfamily = "sans") +
tm_shape(tol.cities) + tm_symbols(shape = 2, col = "red", size = 0.20) +
tm_credits("Data source: UPRA (2020)", fontface = "bold") +
tm_layout(main.title = "Produccion de oleaginosas en 2020",
main.title.fontface = "bold.italic",
legend.title.fontfamily = "monospace") +
tm_scale_bar(position = c("left", "bottom"))
tmap_mode("view")## tmap mode set to interactive viewingcereal_map## Credits not supported in view mode.## Symbol shapes other than circles or icons are not supported in view mode.