Ejercicio: Mapa de aptitud para el cultivo de caña

Actividad

  • Utilizando los datos de clima de línea base a nivel global, genere un código en R que permita construir mapas de aptitud en términos climáticos para la caña de azúcar (con base en los rangos óptimos). Grafique los mapas con una escala de colores adecuada.

  • Identifique 2 o 3 países con áreas de alto potencial para la caña de azúcar y realice un corte para estas zonas con el shape global. Grafique los mapas con una escala de colores adecuada.

  • Identificar algunos puntos (2 o 3) al azar en la región del valle del cauca (use google maps) y extraer la información de clima. Grafique las series de tiempo de temperatura y precipitación.

  • Por medio de alguna métrica de similaridad (ejemplo: distancia euclidiana) genere un código en R que permita identificar mapas de similaridad a nivel global para los sitios identificados en 3. Grafique los mapas con una escala de colores adecuada. Compare los mapas generados por ambas aproximaciones y concluya.

#Cargar las librerias a necesitar
library(raster)
library(rgdal)
library(terra)
library(sf)
library(dplyr)
library(RColorBrewer)
library(plotly)
library(geosphere)
library(leaflet)
library(rnaturalearthdata)
library(rnaturalearth)

Mapas de aptitud en términos climáticos para la caña de azúcar

Primero se obtienen los datos de precipitaciones mensuales en el mundo, luego los datos de temperaturas mensuales mundiales. Teniendo en cuenta que las condiciones óptimas en cuanto a precipitación y temperatura para el cultivo de caña son:

  • Precipitación: 1500mm a 2500mm anuales

  • Temperatura: 20°c a 35°c

Se generó un mapa global promedio de las condiciones óptimas de temperatura y precipitación para le cultivo de caña de azúcar y doce mapas globales que cumplen con esas dos condiciones por mes.

#mapa de precipitacion mundial mensual
ruta=list.files("C:/Users/dolli/Documents/ESTUDIO/Maestria ciencia de datos/Segundo_semestre/electiva/wc2.1_10m_prec/",full.names = TRUE,pattern = ".tif")
prec_year=stack(ruta)
names(prec_year)=month.name

#mapa de temperatura mundial menusales
ruta1=list.files("C:/Users/dolli/Documents/ESTUDIO/Maestria ciencia de datos/Segundo_semestre/electiva/wc2.1_10m_tavg/",full.names = TRUE,pattern = ".tif")
temp_year=stack(ruta1)
names(temp_year)=month.name

# los mapas con precipitacion optima
prec_year_opt=prec_year>=125 & prec_year<=208

# los mapas con temperatura optima
temp_year_opt=temp_year>=20 & temp_year<=35

# condiciones optimas
aptitud_cana = temp_year_opt & prec_year_opt

#promedio de las condiciones optimas cañas mensuales
prom_aptitud_cana = mean(aptitud_cana)
colors = colorRampPalette(brewer.pal(9, "YlGn"))(100)
plot(prom_aptitud_cana, col = colors, 
     main = "Mapa promedio de aptitud para la caña de azúcar")

# Visualizar los mapas con escala de colores adecuada
plot(aptitud_cana, col = colors)

Como se observa en los mapas globales mensuales y en el mapa promedio de aptitud óptima para la caña de azúcar, se puede concluir que los países ubicados cerca de la línea del ecuador son los más propensos a cumplir con estas condiciones durante gran parte del año. Dado que la caña es un cultivo de ciclo largo, esto resultaría beneficioso para mantener una buena producción.

Mapas de tres países con aptitud en términos climáticos para la caña de azúcar

Primero se obtiene el shape file de los limites administrativos de todos lo países en el mundo, esto se cruza con nuestro mapa promedio de aptitud climática para caña de azúcar, y de allí se obtienen el listado de los países de mayor aptitud, para así elegir a los tres primeros.

world = ne_countries(scale = "medium", returnclass = "sf")

# Convertir world a SpatVector
world_terra = as(world, "SpatVector")

#seleccion de paises 
# convertir el raster a puntos
puntos_aptitud = as.data.frame(rasterToPoints(prom_aptitud_cana), xy=TRUE)

#filtrar puntos con alta aptitud
puntos_alta_aptitud = puntos_aptitud[puntos_aptitud$layer > 0.7, ]

#crear un objeto Spatvector con los puntos de alta aptitud
puntos_vect = vect(puntos_alta_aptitud, geom = c('x','y'), crs = "EPSG:4326")

# Extraer el nombre del país correspondiente para cada punto en 'puntos_vect'
extraer_data =  extract(world_terra, puntos_vect)
puntos_vect$country = extraer_data$admin

# Usar aggregate() para calcular la aptitud media por país
aptitud_media_paises = aggregate(layer ~ country, data = puntos_vect, FUN = mean, na.rm = TRUE)

# Ordenar los países por aptitud media en orden descendente
aptitud_ordenada = aptitud_media_paises[order(-aptitud_media_paises$layer), ]

# Seleccionar los tres países con mayor aptitud
top_3_paises = head(aptitud_ordenada$country, 3)
print(top_3_paises)
## [1] "French Polynesia" "Malaysia"         "Peru"

Debido a que la Polinesia Francesa, es un país conformado por un conjunto de islas, su visualización no se hace tan atractiva por la resolución de pixeles que se manejan, por tal motivo se decidió remplazar ese país por Brasil, que es el cuarto país con mas aptitud promedio para caña.

Después de elegir los tres países, vamos a realizar un corte de nuestro mapa mundial de aptitud para caña de azúcar, utilizando las divisiones administrativas de los países elegidos, esta información viene dentro de los atributos del shape global, obteniendo un mapa promedio de condiciones ópticas de caña para cada país elegido y 12 mapas mensuales de esas condiciones por país.

# Filtrar el polígono de Malaysia
malasia_poligon = subset(world_terra, world_terra$admin == "Malaysia")

# Convertir el SpatVector a sf
malasia_poligon_sf = st_as_sf(malasia_poligon)

## Obtiene el bounding box del mapa de Papua Nueva Guinea
malasia_poligon_bbox = st_bbox(malasia_poligon_sf)

# Crear un objeto Extent usando las coordenadas del bounding box
malasia_poligon_raster = extent(c(malasia_poligon_bbox["xmin"], malasia_poligon_bbox["xmax"],
                                 malasia_poligon_bbox["ymin"], malasia_poligon_bbox["ymax"]))

#recortar
malasia_raster = crop(prom_aptitud_cana, malasia_poligon_raster)
# Masking
malasia_raster = terra::mask(malasia_raster, malasia_poligon_sf)
plot(malasia_raster, main = "Mapa de aptitud promedio de caña de azúcar de Malaysia y por mes")

#recorte mensual
malasia_raster_year = crop(aptitud_cana, malasia_poligon_raster)
# Masking
malasia_raster_year = terra::mask(malasia_raster_year, malasia_poligon_sf)
plot(malasia_raster_year)

# Filtrar el polígono de Peru
peru_poligon <- subset(world_terra, world_terra$admin == "Peru")

# Convertir el SpatVector a sf
peru_poligon_sf = st_as_sf(peru_poligon)

## Obtiene el bounding box del mapa de Papua Nueva Guinea
peru_poligon_bbox = st_bbox(peru_poligon_sf)

# Crear un objeto Extent usando las coordenadas del bounding box
peru_poligon_raster = extent(c(peru_poligon_bbox["xmin"], peru_poligon_bbox["xmax"],
                                peru_poligon_bbox["ymin"], peru_poligon_bbox["ymax"]))

#recortar con mapa promedio
peru_raster = crop(prom_aptitud_cana, peru_poligon_raster)
# Masking
peru_raster = terra::mask(peru_raster, peru_poligon_sf)
plot(peru_raster, main = "Mapa de aptitud promedio de caña de azúcar de Perú y por mes")

#recorte con mapa mensual
peru_raster_year = crop(aptitud_cana, peru_poligon_raster)
# Masking
peru_raster_year = terra::mask(peru_raster_year, peru_poligon_sf)
plot(peru_raster_year)

# Filtrar el polígono de Brazil
brasil_poligon = subset(world_terra, world_terra$admin == "Brazil")

# Convertir el SpatVector a sf
brasil_poligon_sf = st_as_sf(brasil_poligon)

## Obtiene el bounding box del mapa de brasil
brasil_poligon_bbox = st_bbox(brasil_poligon_sf)

# Crear un objeto Extent usando las coordenadas del bounding box
brasil_poligon_raster = extent(c(brasil_poligon_bbox["xmin"], brasil_poligon_bbox["xmax"],
                               brasil_poligon_bbox["ymin"], brasil_poligon_bbox["ymax"]))

#recortar con mapa promedio
brasil_raster = crop(prom_aptitud_cana, brasil_poligon_raster)
# Masking
brasil_raster = terra::mask(brasil_raster, brasil_poligon_sf)
plot(brasil_raster, main = "Mapa de aptitud promedio de caña de azúcar de Brasil y por mes")

#recorte con mapa mensual
brasil_raster_year = crop(aptitud_cana, brasil_poligon_raster)
# Masking
brasil_raster_year = terra::mask(brasil_raster_year, brasil_poligon_sf)
plot(brasil_raster_year)

En los tres países elegidos y representados a partir del mapa global, se observa que en ciertas regiones de su territorio se cumplen durante todo el año las condiciones óptimas para el cultivo de caña de azúcar. Esto se refleja en las cifras agrícolas de Brasil, que se posiciona como el líder indiscutible en la producción de caña, seguido de Perú, que también cuenta con una industria significativa. En contraste, Malasia, a pesar de tener condiciones climáticas favorables, ha optado por priorizar cultivos más rentables y estratégicos, como el aceite de palma, debido a una combinación de factores económicos, estratégicos y de disponibilidad de tierras. En conclusión, los mapas creados confirman las condiciones favorables de estos países para el cultivo de la caña de azúcar

Graficas de series de tiempo de temperatura y precipitación para tres puntos en el Valle del Cauca

Utilizando Google Maps, seleccionamos aleatoriamente tres puntos ubicados en el departamento del Valle del Cauca. Con sus coordenadas, extrajimos la información de temperatura y precipitación correspondiente a cada punto durante los 12 meses del año. Luego, graficamos las series de tiempo de cada ubicación para analizar su comportamiento a lo largo del año y compararlas entre sí, así como con los mínimos de temperatura y precipitación requeridos para el cultivo de caña de azúcar. Este análisis nos permitió observar la posible similitud entre los puntos en relación con las condiciones óptimas para el cultivo.

#primer_punto = 3.294676, -76.58812
#segundo_punto = 3.3347311337105325, -76.26559111245042
#tercer_punto = 3.6836203468003044, -76.40128477914159

##extraer de una coordenada punto 1 de una imagen
loc1=cbind(-76.58,3.29)
prec1 = extract(prec_year,loc1)
temp1 = extract(temp_year, loc1)
aptc1 = extract(aptitud_cana, loc1)

##extraer de una coordenada punto 2 de una imagen
loc2=cbind(-76.26,3.33)
prec2 = extract(prec_year,loc2)
temp2 = extract(temp_year, loc2)
aptc2 = extract(aptitud_cana, loc2)

##extraer de una coordenada punto 3 de una imagen
loc3=cbind(-76.40,3.68)
prec3 = extract(prec_year,loc3)
temp3 = extract(temp_year, loc3)
aptc3 = extract(aptitud_cana, loc3)

meses = month.abb

# Crear el DataFrame
df_clima3 = data.frame(
  latitud = rep(c(loc1[,2], loc2[,2], loc3[,2]), each = 12),
  longitud = rep(c(loc1[,1], loc2[,1], loc3[,1]), each = 12),
  Mes = rep(meses, 3),  # Repetir los nombres de los meses
  Temperatura = c(temp1[1,], temp2[1,], temp3[1,]),
  Precipitacion = c(prec1[1,], prec2[1,], prec3[1,]),
  Aptitud_cana_p = c(aptc1[1,], aptc2[1,], aptc3[1,]),
  Punto = rep(c("Punto 1", "Punto 2", "Punto 3"), each = 12))

# Definir los nombres de los meses en español
meses_espanol = c("Enero", "Febrero", "Marzo", "Abril", "Mayo", 
                   "Junio", "Julio", "Agosto", "Septiembre", 
                   "Octubre", "Noviembre", "Diciembre")

# Actualizar data frame para incluir los meses en español
df_clima3$Mes = meses_espanol

# Convertir la columna 'mes' a un factor con el orden correcto
df_clima3$Mes <- factor(df_clima3$Mes, levels = meses_espanol, ordered = TRUE)

# Graficar las lineas (Temperatura)
fig <- plot_ly(df_clima3, 
               x = ~Mes, 
               y = ~Temperatura,
               color = ~Punto,
               type = 'scatter', 
               mode = 'lines + markers', 
               name = ~Punto) 
fig = fig %>%
  add_trace(y = rep(20,36),  
            x = df_clima3$Mes,
            mode = 'lines',
            name = '20 Grados',
            line = list(color = 'red', dash = 'dash'),
            showlegend = FALSE)  
fig %>%
  layout(title = "Series de Tiempo de Temperatura en el Valle del Cauca",
         xaxis = list(title = "Mes"),
         yaxis = list(title = "Temperatura promedio (°C)"))
# Graficar las lineas (Precipitacion)
fig <- plot_ly(df_clima3, 
               x = ~Mes, 
               y = ~Precipitacion,
               color = ~Punto,
               type = 'scatter', 
               mode = 'lines + markers', 
               name = ~Punto) 
fig = fig %>%
  add_trace(y = rep(128,36),  
            x = df_clima3$Mes,
            mode = 'lines',
            name = '128 mm',
            line = list(color = 'red', dash = 'dash'),
            showlegend = FALSE)  
fig %>%
  layout(title = "Series de Tiempo de Precipitación en el Valle del Cauca",
         xaxis = list(title = "Mes"),
         yaxis = list(title = "Preciitación promedio (mm)"))

Como se observa en las gráficas de series de tiempo de temperatura y precipitación, el punto dos presenta las mejores condiciones para el establecimiento de caña de azúcar. Este punto se aleja de los mínimos de temperatura necesarios para el cultivo, y en ocho meses diferentes, su precipitación se encuentra por encima o muy cerca del mínimo mensual requerido. Esto permite que, en promedio, compense mejor estas dos variables climáticas en comparación con los otros puntos.

Además, se diseñó un mapa interactivo de los tres puntos seleccionados, lo que permite realizar comparaciones mensuales de temperatura y precipitación para cada uno de ellos.

### mapa interativo para el punto 3

# Definir paletas de colores para temperatura y precipitación
paleta_temp <- colorNumeric(palette = "viridis", domain = df_clima3$Temperatura)
paleta_prec <- colorNumeric(palette = "plasma", domain = df_clima3$Precipitacion)

# Crear un mapa interactivo temperatura
map_temp = leaflet(df_clima3) %>%
  setView(lng = mean(df_clima3$longitud), lat = mean(df_clima3$latitud), zoom = 6) %>%
  addTiles()  # Añadir el fondo de mapa base

# Título del mapa de temperatura
map_temp <- map_temp %>%
  addControl("Mapa de Temperatura", position = "topright")  

# Añadir marcadores para cada punto
map_temp = map_temp %>%
  addCircleMarkers(~longitud, ~latitud, 
                   color = ~paleta_temp(Temperatura), 
                   radius = 5,
                   label = ~paste("Punto:", Punto, "<br>Mes:", Mes, "<br>Temperatura:", Temperatura),
                   clusterOptions = markerClusterOptions())

# Agregar líneas de tiempo por punto
for (punto in unique(df_clima3$Punto)) {
  temp_data = df_clima3 %>% filter(Punto == punto)
  
  # Agregar la línea de temperatura
  map_temp =  map_temp %>%
    addPolylines(data = temp_data,
                 lng = ~longitud,
                 lat = ~latitud,
                 color = "red", 
                 weight = 2,
                 group = "Temperatura")}



# Crear un mapa interactivo
map_prec = leaflet(df_clima3) %>%
  setView(lng = mean(df_clima3$longitud), lat = mean(df_clima3$latitud), zoom = 6) %>%
  addTiles()  # Añadir el fondo de mapa base

# Título del mapa de precipitación
map_prec <- map_prec %>%
  addControl("Mapa de Precipitación", position = "topright")  

# Añadir marcadores para cada punto
map_prec = map_prec %>%
  addCircleMarkers(~longitud, ~latitud, 
                   color = ~paleta_prec(Precipitacion), 
                   radius = 5,
                   label = ~paste("Punto:", Punto, "<br>Mes:", Mes, "<br>Precipitación:", Precipitacion),
                   clusterOptions = markerClusterOptions())

# Agregar líneas de tiempo por punto
for (punto in unique(df_clima3$Punto)) {
  prec_data = df_clima3 %>% filter(Punto == punto)
  
  # Agregar la línea de precipitacion
  map_prec =  map_prec %>%
    addPolylines(data = prec_data,
                 lng = ~longitud,
                 lat = ~latitud,
                 color = "blue", 
                 weight = 2,
                 group = "Precipitación")}

# Sincronizar los dos mapas
mapjoin <- leafsync::sync(map_temp, map_prec)

# Mostrar los mapas sincronizados
mapjoin

Mapa con métrica de similaridad (distancia euclidiana) entre los puntos geográficos del Valle del Cauca

Lo primero que se realizó fue normalizar los datos a partir del DataFrame de los puntos extraídos, con el fin de crear una matriz de distancias euclidianas. Esta matriz nos permite obtener una métrica de similitud entre los puntos, lo que facilita la identificación de aquellos que son más similares y aquellos que no lo son. Posteriormente, ubicamos los puntos en un mapa, incorporando sus valores de similitud, lo que nos brinda una perspectiva adicional para comparar los puntos.

# Convertir la variable Aptitud_cana a numérica (0 y 1)
df_clima3 <- df_clima3 %>%
  mutate(
    Aptitud_cana_num = as.numeric(Aptitud_cana_p))

#Normalizar las variables para evitar escalas diferentes en los cálculos de distancia
df_clima3 = df_clima3 %>%
  mutate(
    Temp_norm = scale(Temperatura),
    Prec_norm = scale(Precipitacion),
    Aptitud_norm = scale(Aptitud_cana_num))

# Calcular la distancia euclidiana entre los puntos
distancia = dist(df_clima3 %>% select(Temp_norm, Prec_norm, Aptitud_norm))

# Convertir la matriz de distancia en una tabla de distancias
distancia_matriz = as.matrix(distancia)
rownames(distancia_matriz) = df_clima3$Punto
colnames(distancia_matriz) = df_clima3$Punto

# Obtener el valor medio de similitud (menor distancia) para cada punto
similaridad_media <- apply(distancia_matriz, 1, function(x) mean(x[x != 0]))

# Unir las coordenadas con el valor de similitud
df_puntos = df_clima3 %>%
  group_by(Punto) %>%
  summarise(Latitud = first(latitud),
            Longitud = first(longitud)) %>%
  mutate(Similitud = similaridad_media[Punto])

# Crear un mapa de leaflet coloreado por similitud
paleta = colorNumeric(palette = "viridis", domain = df_puntos$Similitud)

mapa = leaflet(df_puntos) %>%
  addTiles() %>%
  addCircleMarkers(
    ~Longitud, ~Latitud,
    radius = 10,
    color = ~paleta(Similitud),
    stroke = FALSE,
    fillOpacity = 0.8,
    popup = ~paste0("Punto: ", Punto, "<br>",
                    "Similitud: ", round(Similitud, 2))
  ) %>%
  addLegend(pal = paleta, values = ~Similitud, title = "Similitud")

# Mostrar el mapa
mapa

De acuerdo con lo que se muestra en el mapa, los tres puntos analizados son distintos entre sí. El punto 1, con una similitud de 1.88, y el punto 2, con una similitud de 1.95, son los menos similares, lo que indica que sus condiciones de temperatura y precipitación difieren de manera significativa. En cambio, el punto 3, con una similitud de 1.93, presenta cierta semejanza tanto con el punto 1 como con el punto 2, aunque se inclina más hacia el segundo. Esto sugiere que sus condiciones de temperatura y precipitación son más similares a las del punto 2, pero no idénticas.

  • Conclusiones

Para un análisis detallado y temporal, el mapa de series de tiempo resulta ser más útil, ya que permite observar las variaciones mensuales y cómo estas se relacionan con los requisitos específicos del cultivo de caña de azúcar.

Pero si se necesita una comparación general de similitudes, el mapa de distancias euclidianas facilita la visualización rápida de qué puntos son más similares y cuáles presentan diferencias significativas en sus condiciones climáticas.

Para concluir, ambos mapas son complementarios, puesto que al comparar las series de tiempo de los tres puntos con el mapa de distancias euclidianas, ambos se complementan y evidencian que el punto 2 presenta las condiciones más favorables para el cultivo de caña de azúcar.