💡 Datos para la actividad:

📝 ACTIVIDAD DE EVALUACIÓN CONTINUA

Los capítulos 3, 4, 5 y 6 de la asignatura revisan los principios del Análisis Exploratorio Espacial, combinando contenidos teóricos y prácticos. En relación a estos últimos, los alumnos deberán entregar como evaluación continua una serie de actividades prácticas que, a final de curso, conformarán un pequeño dossier. Este dossier contendrá material gráfico, en forma de mapas y de gráficos estadístico, y numérico, en forma de tablas conteniendo tanto frecuencias como estadísticos calculados por el alumnos.

Además, a final de curso, los contenidos antes citados de este dossier se acompañarán de un análisis y descripción de las peculiaridades territoriales de los dos países elegidos por el alumno (que deberán estar ubicados en diferentes continentes).

PREPARACIÓN DE LOS DATOS

Comenzaremos la actividad estableciendo la carpeta de trabajo.

Los paquetes necesarios para esta actividad son los siguientes

library(rgeoboundaries)
library(geodata)
library(sf)
library(terra)
library(tmap)
library(tidyverse)

DESCARGA DE LA INFORMACIÓN EN FORMATO VECTORIAL

En primer lugar, hay que descargar el objeto vectorial que corresponde a las fronteras de los países elegidos. En este caso, se descargará la información correspondiente a Austria con la función gb_adm0 del paquete rgeoboundaries, y se simplificará eliminando las variables no necesarias (se conserva sólo la variable que ocupa la tercera columna).

austria_sf <- gb_adm0("Austria")
austria_sf <- austria_sf[3]

DESCARGAS DE LA INFORMACIÓN EN FORMATO RÁSTER.

Las capas que serán utilizadas en esta actividad y en las sigientes son:

Capa de usos de suelo.

La capa de usos de suelo utilizada es la Global Land Cover 2000. Está basada en el “FAO Land Cover Classification System”, identificando originalmente 23 categorías con una resolución de 1 km. Para reducir su peso y complejidad, se ha remuestreado a una resolución de 10 km y se han agrupado las clases originales en 10 categorías

Categoría	Denominación	Colores
1	Superficie Forestal	“#006E00”
2	Manglares	“#9B82E6”
3	Otra vegetación natural	“#E664E6”
4	Pastos	“#B4E664”
5	Humedales	“#B4FEF0”
6	Cultivos	“#F06432”
7	Suelo desnudo	“#C8C8C8”
8	Agua	“#A6CEE3”
9	Glaciares	“#FFFFFF”
10	Zonas artificiales	“#FF0000”

En primer lugar, se carga el fichero ráster que continene la clasificación de usos de suelo de todo el planeta, con una resolución de 1 km.

borr_usos <- rast("D:/G2040/TEMA_3_ESDA_Dossier/datos/usos_5m_10clases.tif")

Para extraer la información correspondiente al país que hemos seleccionado, primero se debe recortar la imagen (función crop()) y luego convertir todos los píxeles situados en el exterior en NoData (función mask()). Estas funciones requieren la misma sintaxis: el fichero ráster a recortar y limpiar y el fichero vectorial con los límites del país.

borr_usos <- crop(borr_usos, austria_sf)
at_usos <- mask(borr_usos, austria_sf)

Para conocer cuántos usos diferentes aparecen en nuestro país acudimos a la siguiente función:

unique(at_usos)

##   glc2000_v1_1
## 1            1
## 2            4
## 3            6
## 4            8
## 5            9
## 6           10

Este tipo de variables discretas (categóricas) se representa en tmap añadiendo a la función tm_raster() el argumento style = "cat". Las categorías de usos de suelo presentes en Austria corresponden a

tm_shape(at_usos) + 
  tm_raster(style= "cat", 
            labels = c("Superficice Forestal", "Pastos", "Cultivos", "Agua", "Glaciares", "Zonas artificiales"), 
            palette = c("#006E00", "#B4E664", "#F06432", "#A6CEE3", "#FFFFFF", "#FF0000"), 
            title="Usos del suelo en Austria") +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Usos del suelo", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

Para proceder a los análisis estadísticos posteriores es necesario convertir los objetos cartográficos en dataframes.

at_usos_df <- as.data.frame(at_usos, xy = TRUE, na.rm = TRUE) %>% 
  rename(usos = glc2000_v1_1) %>% mutate(usos = as.factor(round(usos)))

Para no cargar el Global Environment con objetos que no se utilizan, se procede a su eliminación:

rm(borr_usos)

Clasificación climática

La siguiente base de datos incluye una clasificación de los climas del planeta siguiendo el método de [Köppen-Geiger] (https://koeppen-geiger.vu-wien.ac.at/present.htm]). Esta base de datos proporciona la distribución espacial de los diferentes tipos de clima en todo el planeta con una resolución espacial de 0,5 grados, representativa del periodo 1951-2000. Al igual que la anterior, ha sido remuestreada para adaptase a una resolución de 5 minutos y simplificada en un total de 9 grandes conjuntos climáticos.

Las categorías y los colores asociados son los siguientes

Categoría	Denominación	Colores
1	Clima tropical (A)	“#FF0000”
2	Clima desértico (B)	“#CC8D14”
3	Clima templado húmedo (Cf)	“#007800”
4	Clima templado subhúmedo (Cw)	“#96FF00”
5	Clima mediterráneo (Cs)	“#BEBE00”
6	Clima continental húmedo (Df)	“#550055”
7	Clima continental frío (Ds)	“#646464”
8	Clima continental mediterraneizado (Dw)	“#6E28B4”
9	Clima fríos Polar/Alta montaña (E)	“#6496FF”

En primer lugar, se cargará el objeto ráster.

borr_clima <- rast("D:/G2040/TEMA_3_ESDA_Dossier/datos/koppen_5m_9clases.tif")

Como en el caso anterior, debe recortarse y enmascararse el país de elección

borr_clima <- crop(borr_clima, austria_sf)
at_clima <- mask(borr_clima, austria_sf)

Para conocer cuántos climas pueden hallarse en ese país, se acude a la función unique()

unique(at_clima)

##   layer
## 1     3
## 2     6
## 3     9

Los tipos de clima también constituyen una variable discretas (categóricas), como el caso de los usos de suelo. Los tipos de clima que aparecen en Austra son los siguientes.

tm_shape(at_clima) + 
  tm_raster(style= "cat", 
            labels = c("Clima templado húmedo (Cf)", "Continental húmedo (Df)", "Alta montaña (E)"), 
            palette = c("#007800", "#550055", "#6496FF"), 
            title="Tipos de climas") +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Climas", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

Como se realizó en el caso anterior, convertimos el objeto at_clima en un dataframe.

at_clima_df <- as.data.frame(at_clima, xy = TRUE, na.rm = TRUE) %>% 
  rename(clima = layer) %>% mutate(clima = as.factor(round(clima)))

Para no llenar el Global Environment de objetos sin uso, se eliminan los siguientes:

rm(borr_clima)

Clasificación formas de relieve

La última variable cualitativa que se incorpora es una clasificación de [formas de relieve del planeta] (https://esdac.jrc.ec.europa.eu/content/global-landform-classification/). Esta base de datos ha sido elaborada por ESDAC - European Commission.

borr_relieve <- rast("D:/G2040/TEMA_3_ESDA_Dossier/datos/landform_5m_9clases.tif")

En este caso, las categorías y los colores asociados son los siguientes

Categoría	Denominación	Colores
1	Llanuras	“#A7DFD2”
2	Depresiones	“#1F78B4”
3	Mesetas	“#EFEBC0”
4	Zonas onduladas	“#CAB982”
5	Montañas	“#AA8753”

Recorte y enmascaramiento del relieve

borr_relieve <- crop(borr_relieve, austria_sf)
at_relieve <- mask(borr_relieve, austria_sf)

Para conocer cuántos usos diferentes aparecen en nuestro país acudimos a la siguiente función:

unique(at_relieve)

##   landform_esdac
## 1              1
## 2              2
## 3              3
## 4              4
## 5              5

Para representar el mapa de formas de relieve de Austria:

tm_shape(at_relieve) + 
  tm_raster(style= "cat", 
            labels = c("Llanuras", "Depresiones", "Mesetas", "Zonas onduladas", "Montañas"), 
            palette = c("#A7DFD2", "#1F78B4", "#EFEBC0", "#CAB982", "#AA8753"), 
            title="Formas de relieve") +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, tipos de relieve", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

Como se realizó en el caso anterior, convertimos el objeto at_relieve en un dataframe.

at_relieve_df <- as.data.frame(at_relieve, xy = TRUE, na.rm = TRUE) %>% 
  rename(relieve = landform_esdac) %>% mutate(relieve = as.factor(round(relieve)))

Por último se elimina lo que no es necesario.

rm(borr_relieve)

Modelo digital de elevaciones.

Este modelo será extraido de la base de datos contenida en el paquete geodata. Se comenzará cargando el modelo global con la resolución (5 minutos) adecuada.

borr_mdt <- elevation_global(res = 5, path=tempdir())

Posteriormente, se recorta y enmascara la zona de interés (Austria, en este caso)

borr_mdt <- crop(borr_mdt, austria_sf)
at_altitud <- mask(borr_mdt, austria_sf)

Un mapa de isolíneas es la representación en dos dimensiones de una superficie estadística suavizada, a través de elementos lineales, denominados isolíneas. Las isolíneas, como son las curvas de nivel o las isobaras, etc…), son líneas que unen puntos en donde una variable toma un mismo valor. Un mapa topográfico requiere un fichero original en formato raster que incluye una variable continua (en este caso, la altitud); tm_raster(). El argumento style = cont es más adecuado para representar fenómenos que varían de manera progresiva sobre el territorio.

tm_shape(at_altitud) +
   tm_raster(title = "Altitud (mnm):", 
             style = "cont", 
             palette = "-Spectral") +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Hipsometría", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

La creación de un mapa de isopletas requiere crearlas primero con la función tm_iso() y posteriormente cartografiarlas. La creación de isopletas se realiza con la función as.contour del paquete terra, siendo luego transformadas al formato sf. .

isopletas <- as.contour(at_altitud)
isopletas <- sf::st_as_sf(isopletas)

La representación cartográfica se realiza mediante la función tm_iso()

tm_shape(at_altitud) +
   tm_raster(title = "Altitud (mnm):", 
             style = "cont", 
             palette = terrain.colors(9)) +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_shape(isopletas) +
  tm_iso() + 
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Hipsometría", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

El mapa anterior será convertido en dataframe para su posterior análisis estadístico.

at_altitud_df <- as.data.frame(at_altitud, xy = TRUE, na.rm = TRUE) %>% 
  rename(altitud = wc2.1_5m_elev)

Por último, se eliminan los materiales no necesarios.

rm(borr_mdt)

Datos climáticos globales.

La base de datos WorldClim proporciona datos de diferentes variables climáticas con varias resoluciones espaciales. Se puede extraer directamente desde internet con el paquete **geodata*. A continuación se extraerá la precipitación media anual.

borr_pp <- worldclim_global(var="prec", res = 5, path=tempdir())

Recorte y enmascaramiento del objeto gráfico

at_pp <- crop(borr_pp, austria_sf)
at_pp <- mask(at_pp, austria_sf)
at_pp

## class       : SpatRaster 
## dimensions  : 32, 92, 12  (nrow, ncol, nlyr)
## resolution  : 0.08333333, 0.08333333  (x, y)
## extent      : 9.5, 17.16667, 46.33333, 49  (xmin, xmax, ymin, ymax)
## coord. ref. : lon/lat WGS 84 (EPSG:4326) 
## source(s)   : memory
## names       : wc2.1~ec_01, wc2.1~ec_02, wc2.1~ec_03, wc2.1~ec_04, wc2.1~ec_05, wc2.1~ec_06, ... 
## min values  :          19,          22,          28,          37,          53,          64, ... 
## max values  :         117,         106,         135,         160,         162,         195, ...

El fichero original contiene varias capas, cada una de ellas correspondiente a un mes del año. Para calcular el total anual de precipitaciones debe aplicarse la función app() del paquete terra, a través de la cual se pueden calcular diversos estadísticos.

at_pp <- terra::app(at_pp, sum)

A continuación se representa el mapa de la precipitación anual.

tm_shape(at_pp) +
   tm_raster(title = "mm", 
             style = "cont", 
             palette = "Spectral") +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Precipitación anual", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

Como se realizó en el caso anterior, convertimos el objeto at_relieve en un dataframe.

at_pp_df <- as.data.frame(at_pp, xy = TRUE, na.rm = TRUE) %>% 
  rename(precip_anual = sum)

Por último, se borra lo que no es necesario.

rm(borr_pp)

Población mundial.

Esta base de datos incluye la población total que habita dentro de cada uno de los píxeles. La fuente es la Gridded Population of the World (GPW).

En primer lugar, se descarga la información.

borr_pob <- rast("D:/G2040/TEMA_3_ESDA_Dossier/datos/densidad_poblacion_5m.tif")

Además, se recorta y enmascara el área de trabajo.

at_pob <- crop(borr_pob, austria_sf)
at_pob <- mask(at_pob, austria_sf)

A continuación se representa el mapa con la población de Austria.

tm_shape(at_pob) +
   tm_raster(title = "Habitantes km2", 
             style = "cont", 
             palette = c("lightyellow", "orange", "brown")) +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Población", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

Como se realizó en el caso anterior, convertimos el objeto at_pob en un dataframe.

at_pob_df <- as.data.frame(at_pob, xy = TRUE, na.rm = TRUE)

Eliminación de material no necesario.

rm(borr_pob)

Huella ecológica

La huella ecológica es una estimación de la presión humana, tanto directa como indirecta, en el medio natural. Es un indicador complejo obtenido a partir de la combinación de diferents variables, como las vías de comunicación, la densidad de la población, infraestructuras eléctricas, superfiice cultivada etc. Se expresa en uan escala de 0 (huella reducida) a 50 (huella alta). Primero se cargará el fichero.

borr_huella <- rast("D:/G2040/TEMA_3_ESDA_Dossier/datos/huella_ecologica_5m.tif")

Como en otras ocasiones, se recorta y enmascara.

at_huella <- crop(borr_huella, austria_sf)
at_huella <- mask(at_huella, austria_sf)

A continuación se representa en un mapa.

tm_shape(at_huella) +
   tm_raster(title = "Valores", 
             style = "cont", 
             palette= "Reds") +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Huella ecológica", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

Como se realizó en el caso anterior, convertimos el objeto at_relieve en un dataframe.

at_huella_df <- as.data.frame(at_huella, xy = TRUE, na.rm = TRUE)

Eliminación de materiales no necesarios.

rm(borr_huella)

Accesibilidad a un gran puerto

El precio de exportación de la producción agraria e industrial de un país depende en gran medida de su fácil acceso a los grandes puertos, desde donde se embarcan los productos. A continuación, se descargará un fichero en formato ráster con el tiempo de viaje a un gran puerto.

borr_acces <- rast("D:/G2040/TEMA_3_ESDA_Dossier/datos/accesibilidad_puerto_5m.tif")

Recortamos y enmascaramos

at_acces <- crop(borr_acces, austria_sf)
at_acces <- mask(at_acces, austria_sf)

Se representa gráficamente

tm_shape(at_acces) +
   tm_raster(title = "Distancia (km) al puerto más cercano", 
             style = "cont", 
             palette= c("darkgreen", "white", "red")) +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Accesibilidad marítima", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

Como se realizó en el caso anterior, convertimos el objeto at_relieve en un dataframe.

at_acces_df <- as.data.frame(at_acces, xy = TRUE, na.rm = TRUE)

Eliminación de los materiales no necesarios.

rm(borr_acces)

Superficie cultivada

El último parámetro que será estudiado es la superficie cultivada. A continuación, se descargará un fichero en formato ráster.

borr_cult <- rast("D:/G2040/TEMA_3_ESDA_Dossier/datos/superficie_cultivada_5m.tif")

Recortamos y enmascaramos

at_cult <- crop(borr_cult, austria_sf)
at_cult <- mask(at_cult, austria_sf)

Se representa gráficamente

tm_shape(at_cult) +
   tm_raster(title = "% superficie cultivada", 
             style = "cont", 
             palette= c("darkgreen", "white", "red")) +
tm_shape(austria_sf) +
  tm_borders(col = "black", 
             lwd = 1.5) +
tm_layout(legend.outside = TRUE, 
          frame = FALSE, 
          main.title = "Austria, Superficie cultivada", 
          legend.position = c("left", "top")) +
tm_compass(position = c("right", "bottom")) +
tm_scale_bar(position = c("left", "bottom"), text.size = 0.7)

Como se realizó en el caso anterior, convertimos el objeto at_relieve en un dataframe.

at_cult_df <- as.data.frame(at_cult, xy = TRUE, na.rm = TRUE)

Eliminación de los materiales no necesarios:

rm(borr_cult)

CREACIÓN DE UN ÚNICO DATAFRAME CON TODA LA INFORMACIÓN

at_df <- merge(at_usos_df, at_clima_df, by = c("x", "y"))
at_df <- merge(at_df, at_relieve_df, by = c("x", "y"))  
at_df <- merge(at_df, at_altitud_df, by = c("x", "y"))
at_df <- merge(at_df, at_huella_df, by = c("x", "y"))
at_df <- merge(at_df, at_pp_df, by = c("x", "y"))
at_df <- merge(at_df, at_acces_df, by = c("x", "y"))
at_df <- merge(at_df, at_cult_df, by = c("x", "y"))
at_df <- merge(at_df, at_pob_df, by = c("x", "y"))

GRABACIÓN DE LA INFORMACIÓN

Los resultados de esta sesión se grabarán en dos ficheros.

Por un lado, se grabará todo el material cartográfico en un único fichero *.RData por si fuera necesario un uso posterior.

save(at_sf, at_clima, at_cult, at_huella, at_pob, at_pp, at_relieve, at_usos, 
     file = "at_mapas.RData")

Por otro lado, se grabará el dataframe con toda la información en otro fichero.

save(at_df, 
     file = "at_df.RData")

TEMA 2: ESDA: DOSSIER

Domingo Rasilla y Francisco Conde

14 octubre, 2024

PREPARACIÓN DE LOS DATOS

DESCARGA DE LA INFORMACIÓN EN FORMATO VECTORIAL

DESCARGAS DE LA INFORMACIÓN EN FORMATO RÁSTER.

Huella ecológica

Accesibilidad a un gran puerto

Superficie cultivada

CREACIÓN DE UN ÚNICO DATAFRAME CON TODA LA INFORMACIÓN

GRABACIÓN DE LA INFORMACIÓN