CASO 1

### Librerias
library(sf)        # Para manejar shapefiles
library(dplyr)     # Para análisis de datos
library(ggplot2)   # Para visualización de datos
library(sp)
library(geoR)
library(terra)
library(raster)

Objetivos del Informe

Análisis Exploratorio

En esta parte se debe realizar un analisis exploratorio de las variables de estudio Temperatura y Humedad Relativa en el cual se logre identificar el comportamiento (tendencia central, varibilidad). Usar indicadores y graficos.

Carga de datos/información gráfica

area <- "C:/Users/veron/OneDrive/Documents/ANALISIS_GEO/zona/area.shp"
geoR = "C:/Users/veron/OneDrive/Documents/ANALISIS_GEO/datos_geo.RData"
# Cargar los datos geográficos
load(geoR) # Esto cargará las variables contenidas en el archivo .RData
# Cambiamos los nombres usando rename
datos <- datos_geo %>%
  rename(
    Temperatura = Temperature,
    Humedad_Relativa = Relative_Humidity)

Resumen estadístico

# Resumen estadístico de los datos
summary(datos[,1:2])
##   Temperatura    Humedad_Relativa
##  Min.   :16.00   Min.   :28.50   
##  1st Qu.:19.23   1st Qu.:41.65   
##  Median :21.30   Median :63.25   
##  Mean   :22.62   Mean   :58.09   
##  3rd Qu.:26.48   3rd Qu.:68.97   
##  Max.   :34.60   Max.   :91.40

Gráficos (Histogramas)

## Histograma de temperatura

ggplot(datos, aes(x = Temperatura)) + 
  geom_histogram(bins = 30, fill = "red", alpha = 0.7) + 
  theme_minimal() + 
  labs(title = "Histograma de Temperatura", x = "Temperatura", y = "Frecuencia")

En este gráfico se puede observar que la frecuencia de los datos desciende, hay mayor concentración de datos en los valores más pequeños, sobre todo entre los 15° y 20° de temperatura.

## Histograma de humedad relativa

ggplot(datos, aes(x = Humedad_Relativa)) + 
  geom_histogram(bins = 30, fill = "pink", alpha = 0.7) + 
  theme_minimal() + 
  labs(title = "Histograma de Humedad Relativa", x = "Humedad Relativa", y = "Frecuencia")

En este gráfico se observa que los valores están más dispersos entre ellos, hay una concentración mayor de temperatura en valores como 40 y entre 60 y 80

# Explorar los datos cargados
summary(geoR)
##    Length     Class      Mode 
##         1 character character
str(geoR)
##  chr "C:/Users/veron/OneDrive/Documents/ANALISIS_GEO/datos_geo.RData"

Análisis Exploratorio Espacial

En este punto se deben realizar mapas que permitan identificar el comportamiento espacial de las variables (Figura 1.) interpretar el comportamiento observado.

Del anteior gráfico se puede observar que la zona Noroeste presenta mayores valores de temperatura, mientras que la zona sur presenta menores valores de temperatura, además los datos se dispersan de manera uniforme y homogenea.

Análisis descriptivo de la capa shapefile

## [1] POLYGON
## 18 Levels: GEOMETRY POINT LINESTRING POLYGON MULTIPOINT ... TRIANGLE
##       xmin       ymin       xmax       ymax 
## -76.614175   2.380766 -76.612507   2.382726
## Coordinate Reference System:
##   User input: +proj=longlat +datum=WGS84 +no_defs 
##   wkt:
## GEOGCRS["unknown",
##     DATUM["World Geodetic System 1984",
##         ELLIPSOID["WGS 84",6378137,298.257223563,
##             LENGTHUNIT["metre",1]],
##         ID["EPSG",6326]],
##     PRIMEM["Greenwich",0,
##         ANGLEUNIT["degree",0.0174532925199433],
##         ID["EPSG",8901]],
##     CS[ellipsoidal,2],
##         AXIS["longitude",east,
##             ORDER[1],
##             ANGLEUNIT["degree",0.0174532925199433,
##                 ID["EPSG",9122]]],
##         AXIS["latitude",north,
##             ORDER[2],
##             ANGLEUNIT["degree",0.0174532925199433,
##                 ID["EPSG",9122]]]]

En base a lo anterior nos podemos dar cuenta de que la figura constituida es un poligono y tiene 18 niveles difeentes de geometría.

Análisis Exploratorio Espacial

En este punto se deben realizar mapas que permitan identificar el comportamiento espacial de las variables (Figura 1.) interpretar el comportamiento observado.

g1 = ggplot() + 
  geom_polygon(data=zona, aes(x=lat, y=long), color="white", fill="grey80") +
  geom_point(data=datos_geo, aes(x=Latitude, y=Longitude, size=Temperature, color=Temperature)) +
  theme_bw() +
  scale_color_gradient(low="blue3", high="red")
## Regions defined for each Polygons

Caso 2: Estudio Violencia en la Ciudad de Cali

Se proporciona información de los homicidios registrados en un periodo de tiempo en la ciudad de cali y se desea con base en estos datos realizar una estrategia que permita por medio de la identidicación de caracteristicas relevantes a nivel de individuo y espacial cuales son las zonas de mayor prioridad para mejorar la seguridad en la ciudad y las caracteristicas de los eventos.

Usted tiene libertad para realizar este trabajo solo se requiere que en alguno de los puntos realice uso de patrones puntuales.

Datos que tengo:

shape points homicidios shape de comunas densidad por comuna barrio y estratos estaciones de policia y CAIS

Lo primero que se quiere ver es la relación de comunas y homicidios, es decir, poder visualizar si los homicidios se acumulan en algunas comunas en específico. Para hacer esto es necesario utilizar patrones puntuales, y algo a tener en cueta en patrones puntuales es que pueden ser regulares, aleatorios o agregado. Para identificar el tipo de patrón se va a utilizar la función K.

Función K:

#Función K para ver qué tipo de distribución tiene:
#Creación de la ventana de datos
bbox <- st_bbox(homicidios)
ventana <- owin(xrange = c(bbox["xmin"], bbox["xmax"]),
                yrange = c(bbox["ymin"], bbox["ymax"]))
#conversión de datos shp a ppp
#homicidios_ppp <- as.ppp(st_coordinates(homicidios), W = ventana)
homicidios_ppp <- as.ppp(st_coordinates(homicidios), W = ventana, dimyx = c(500, 500))
## Warning: data contain duplicated points
plot(Kest(homicidios_ppp, correction = c("iso"), ratio = TRUE))

#R/ tienen un comportamiento: agregado. Kiso está por arriva de Kpois

En la gráfica anterior Kpois representa un comportamiento aleatorio, a partir de esto se analiza la posición y comportamiento de Kiso; y se puede observar que Kiso está arriba de Kpois, lo que significa que tiene un comportamiento agregado. Esto significa, en el contexto del ejercicio, que los homicidios están agregados en algunas zonas específicas de la ciudad.

Mapa de homicidios por comuna:

A continuación se presenta un mapa donde se puede visualizar el número de homicidios por zonas, de manera acumulada, y también, si se hace zoom se puede ver la jornada (mañana, tarde, noche) en la cual fue cometido cada homicidio.

#Visualización de homicidios agrupados
mapa1_h <- leaflet(homicidios) %>% addTiles() %>% addCircleMarkers(clusterOptions = markerClusterOptions(), label = ~jornada) %>% 
  addControl(html = "<h3>Mapa de Homicidios<h3>", position = "topleft")
mapa1_h <- mapa1_h %>% addPolygons(data = comunas, label =  comunas$COMUNA, fillColor = "transparent", color = "#2258c7", weight = 1, opacity = 1)
mapa1_h

Con el mapa anterior se puede observar que la mayoría de homicidios se acumulan en la zona oriente y en la zona de la ladera oeste. Sin embargo, para saber con exactitud cuáles son las comunas con mayores homicidios,se puede visualizar esta información a través de tablas, donde se muestran las 5 comunas con mayor cantidad de homicidios, también, se muestra el porcentaje de homicidios en la comuna respecto al total de homicidios.

Tabla de comunas con mayores homicidios

#Contar número de homicidios por comuna:
#Se utiliza un join entre homicidios y comunas, para contar los puntos que se encuentran dentro del polígono de cada comuna
join_homi_comu <- st_join(homicidios, comunas, join = st_within)

#Se agrupan por el número de comunas, se cuentan los homicidios y se organizan de mayor a menor
cant_homicidios <- join_homi_comu %>%
  group_by(COMUNA) %>%
  summarise(cant_homicidios = n())

total_homicidios <- sum(cant_homicidios$cant_homicidios)

cant_homicidios <- cant_homicidios %>%
  mutate(porcentaje_homicidios = (cant_homicidios / total_homicidios) * 100)

tabla_cant_homicidios <- cant_homicidios %>%
  st_drop_geometry() %>%
  arrange(desc(cant_homicidios)) %>%
  slice_head(n = 5)
Top 5 Comunas con Mayor Número de Homicidios
Comuna Cantidad de Homicidios Porcentaje de Homicidios
15 552 11.66 %
14 504 10.64 %
13 488 10.30 %
21 398 8.40 %
20 355 7.50 %

Con la esta tabla se puede ver que las comunas con más homicidios son la 15, 14, 13, 21, 20. La suma de los homicidios en estas comunas acumulan más del 40% de los homicidios que ocurren en toda la ciudad, lo que significa que son zonas críticas y de alta prioridad.

Caracterización homicidios

Homicidios por barrios:

A continaución se presentan los homicidios distribuidos por barrios en las 5 comunas identificadas:

#Visualización de homicidios agrupados por barrios
mapa2_h <- leaflet(homicidios) %>% addTiles() %>% addCircleMarkers(clusterOptions = markerClusterOptions(), label = ~jornada) %>% 
  addControl(html = "<h3>Homicidios en barrios<h3>", position = "topleft")
# polígonos de las comunas
mapa2_h <- mapa2_h %>%
  addPolygons(data = comunas_t, fillColor = "orange", color = "white", weight = 2, opacity = 1)
#polígonos de los barrios
mapa2_h <- mapa2_h %>%
  addPolygons(data = barrios_comunas_5, label = ~BARRIO, fillColor = "transparent", color = "#2258c7", weight = 1, opacity = 1)
mapa2_h

Estaciones de policía

Ahora vamos a analizar las estaciones de policía a lo largo de la ciudad. Esta variable es tipo punto entonces empezamos viendo si tiene algún tipo de distribución específica o patrón puntual.

estaciones <- st_transform(estaciones, crs = 4326)
estaciones <- estaciones %>% filter(tipo %in% c("CAI", "ESTPOL"))
#Función K para ver qué tipo de distribución tiene:
#Creación de la ventana de datos
bbox <- st_bbox(estaciones)
ventana <- owin(xrange = c(bbox["xmin"], bbox["xmax"]),
                yrange = c(bbox["ymin"], bbox["ymax"]))
#conversión de datos shp a ppp
#homicidios_ppp <- as.ppp(st_coordinates(homicidios), W = ventana)
estaciones_ppp <- as.ppp(st_coordinates(estaciones), W = ventana, dimyx = c(500, 500))
plot(Kest(estaciones_ppp, correction = c("iso"), ratio = TRUE))

Esta gráfica, a diferencia de los homicidios no presenta un comportamiento altamente agregado, esto se puede decir ya que la gráfica de las estaciones (Kiso) se acerca bastante a la gráfica de aleatoriedad (Kpois). A continuación se observan las estaciones de policiía a lo largo de Cali por comunas

#Visualización de homicidios agrupados
mapa3_h <- leaflet(estaciones) %>% addTiles() %>% addCircleMarkers(clusterOptions = markerClusterOptions()) %>% 
  addControl(html = "<h3>Mapa de Homicidios<h3>", position = "topleft")
mapa3_h <- mapa3_h %>% addPolygons(data = comunas, label =  comunas$COMUNA, fillColor = "transparent", color = "#2258c7", weight = 1, opacity = 1)
mapa3_h
join_esta_comu <- st_join(estaciones, comunas, join = st_within)

#Se agrupan por el número de comunas, se cuentan los homicidios y se organizan de mayor a menor
cant_estaciones <- join_esta_comu %>%
  group_by(COMUNA) %>%
  summarise(cant_estaciones = n())

total_estaciones <- sum(cant_estaciones$cant_estaciones)

tabla_cant_estaciones <- cant_estaciones %>%
  st_drop_geometry() %>%
  arrange(desc(cant_estaciones)) %>%
  slice_head(n = 5)
tabla_cant_estaciones
## # A tibble: 5 × 2
##   COMUNA cant_estaciones
##    <dbl>           <int>
## 1      2               7
## 2     19               4
## 3      3               3
## 4      5               3
## 5      6               3

Con base a lo anterior las entidades de seguridad y policía pueden tomar la decisión de poner más puestos y estaciones en las comunas con más homicidios, es decir, 15, 14, 13, 20, 21. Estas comunas son de mayor prioridad y no tienen la seguridad necearia

Variables de los homicidios:

Sexo:

# Calcular los porcentajes de homicidios por sexo
homicidios_porcentaje_sexo <- homicidios %>%
  group_by(SEXO) %>%
  summarise(cantidad = n(), .groups = "drop") %>%
  mutate(porcentaje = (cantidad / sum(cantidad)) * 100)

# Crear el gráfico de barras
ggplot(homicidios_porcentaje_sexo, aes(x = SEXO, y = porcentaje, fill = SEXO)) +
  geom_bar(stat = "identity") +
  labs(title = "Porcentaje de Homicidios por Sexo",
       x = "Sexo",
       y = "Porcentaje (%)",
       fill = "Sexo") +
  theme_minimal() +
  geom_text(aes(label = paste0(round(porcentaje, 1), "%")), vjust = -0.5) +
  scale_y_continuous(labels = scales::percent_format(scale = 1))

La mayoría de los homicidios se cometen hacia los hombres.

# Contar homicidios por sexo en todas las comunas
# Realizar la unión espacial para asignar cada homicidio a una comuna
homicidios_comunas <- st_join(homicidios, comunas, join = st_within)
homicidios_por_sexo <- homicidios_comunas %>%
  group_by(COMUNA, SEXO) %>%
  summarise(cantidad = n()) %>%
  ungroup()
## `summarise()` has grouped output by 'COMUNA'. You can override using the
## `.groups` argument.
# Gráfico de barras
ggplot(homicidios_por_sexo, aes(x = factor(COMUNA), y = cantidad, fill = SEXO)) +
  geom_bar(stat = "identity", position = "dodge") +
  labs(title = "Cantidad de Homicidios por Sexo en Comunas",
       x = "Comuna",
       y = "Cantidad de Homicidios",
       fill = "Sexo") +
  theme_minimal()

Jornada:

# Calcular los porcentajes de homicidios por sexo
homicidios_porcentaje_jornada <- homicidios %>%
  group_by(jornada) %>%
  summarise(cantidad = n(), .groups = "drop") %>%
  mutate(porcentaje = (cantidad / sum(cantidad)) * 100)

# Crear el gráfico de barras
ggplot(homicidios_porcentaje_jornada, aes(x = jornada, y = porcentaje, fill = jornada)) +
  geom_bar(stat = "identity") +
  labs(title = "Porcentaje de Homicidios por jornada",
       x = "Jornada",
       y = "Porcentaje (%)",
       fill = "Jornada") +
  theme_minimal() +
  geom_text(aes(label = paste0(round(porcentaje, 1), "%")), vjust = -0.5) +
  scale_y_continuous(labels = scales::percent_format(scale = 1))

La jornada en la que más se comenten homicidios es la noche, con más del 40%. Con esta información se peude saber a qué horas de deberían ahcer rondas, así como también pueden influir otras variables como por ejemplo la iluminación de la ciudad.

comunas_seleccionadas <- comunas %>%
  filter(COMUNA %in% c(13, 14, 15, 20, 21))
#join espacial para asignar la comuna a cada homicidio
homicidios_comunas <- st_join(homicidios, comunas_seleccionadas, join = st_within)
#Filtrar homicidios
homicidios_filtrados <- homicidios_comunas %>%
  filter(!is.na(COMUNA))
#cantidad de homicidios por comuna y jornada
homicidios_comuna_jornada <- homicidios_filtrados %>%
  group_by(COMUNA, jornada) %>%
  summarise(cantidad = n(), .groups = "drop")
#porcentaje de homicidios dentro de cada comuna
homicidios_comuna_jornada <- homicidios_comuna_jornada %>%
  group_by(COMUNA) %>%
  mutate(total_comuna = sum(cantidad),
         porcentaje = (cantidad / total_comuna) * 100) %>%
  ungroup()
# 6. Crear el gráfico de barras por comuna
grafico <- ggplot(homicidios_comuna_jornada, aes(x = jornada, y = cantidad, fill = jornada)) +
  geom_bar(stat = "identity", color = "black") +
  geom_text(aes(label = paste0(round(porcentaje, 1), "%")), vjust = -0.5) +
  labs(
    title = "Homicidios por Jornada en Comunas Seleccionadas",
    x = "Jornada",
    y = "Cantidad de Homicidios",
    fill = "Jornada"
  ) +
  facet_wrap(~ COMUNA, scales = "free_y") +  # Un gráfico por comuna
  theme_minimal()

# Mostrar el gráfico
print(grafico)

Días de la semana:

# Calcular los porcentajes de homicidios por sexo
homicidios_porcentaje_dias <- homicidios %>%
  group_by(diasem) %>%
  summarise(cantidad = n(), .groups = "drop") %>%
  mutate(porcentaje = (cantidad / sum(cantidad)) * 100)

# Crear el gráfico de barras
ggplot(homicidios_porcentaje_dias, aes(x = diasem, y = porcentaje, fill = diasem)) +
  geom_bar(stat = "identity") +
  labs(title = "Porcentaje de Homicidios por día de la semana",
       x = "Día se la semana",
       y = "Porcentaje (%)",
       fill = "Día") +
  theme_minimal() +
  geom_text(aes(label = paste0(round(porcentaje, 1), "%")), vjust = -0.5) +
  scale_y_continuous(labels = scales::percent_format(scale = 1))

El día con mas homicidios es el domingo. Nuevamente, se puede analizar cuáles patrones ocurren los domingos, por ejemplo, las personas están en fiestas y tomando bebidas embriagantes, son días para salir de casa, entre otros.

Caso 3.

Caso Identificación de zonas para la caña de azucar Utilizando los datos de clima de línea base a nivel global del siguiente enlace https://www.worldclim.org/data/worldclim21.html, 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. ver video:https://www.youtube.com/watch?v=XANenU2XDQ4&ab_channel=CentroMagis%5BJaverianaCali%5D

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 usando leaflet o mapview de forma interactiva.

CARGA DE INFORMACIÓN

Datos a tener en cuenta (según el video):

Temperatura media: 22.5°c y 28°c Precipitación anual: 1500mm - 3500mm Precipitación mensual: 125mm - 290mm

Carga de información (librerias, datos)

library(raster)
library(sf)
library(leaflet)
library(terra) 
library(mapview)
#Rutas
#Ruta carpeta de datos de temperatura
rutaTem = "C:/Users/veron/Downloads/wc2.1_2.5m_tavg/"

#Ruta carpeta de datos de precipitación
rutaPre = "C:/Users/veron/Downloads/wc2.1_2.5m_prec/"

#Ruta shp global
ruta_paises = "C:/Users/veron/OneDrive/Documents/ANALISIS_GEO/shape_global/g2008_0.shp"

paises <- st_read(ruta_paises)
## Reading layer `g2008_0' from data source 
##   `C:\Users\veron\OneDrive\Documents\ANALISIS_GEO\shape_global\g2008_0.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 534 features and 11 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: -180 ymin: -55.72333 xmax: 180 ymax: 83.62742
## Geodetic CRS:  WGS 84

PREPROCESAMIENTO

El preprocesamiento de los datos consta de ordenar los datos que se tienen para su posterior procesamieto, para el caso 3 el preprocesamiento consta de operaciones sencillas como la suma y el promedio de un conjunto de datos, para poder utilizar la información debidamente.

Preprocesamiento temperatura media:

Una de las condiciones para tener una zona optima para el cultivo de caña de azucar es la temperatura media, para ello, luego de descargar la información en el portal de World clim se realiza la operación de promedio entre los meses de un año, ya que la información que suministra este sitio web está dividido en meses, por lo que se realiza el promedio con todos los meses del año como se observa a continuación:

#TEMPERATURA MEDIA (ANUAL)
# Ruta base de los rasters mensuales
#rutaTem = "C:/Users/hilla/OneDrive/Escritorio/Semestre 5/Analisis de datos espaciales/Taller1_analisis/wc2.1_30s_tavg/"

# Crear una lista con los nombres de los archivos de los 12 meses
archivos_meses = list.files(path = rutaTem, pattern = "wc2.1_2.5m_tavg_.*\\.tif$", full.names = TRUE)

# Cargar todos los rasters en una lista y combinarlos
rasters_mensuales = rast(archivos_meses)

# Calcular el promedio anual de temperatura
prom_Temp = mean(rasters_mensuales)
## |---------|---------|---------|---------|=========================================                                          
# Mostrar el raster del promedio
plot(prom_Temp, main = "Promedio Anual de Temperatura")

Preprocesamiento precipitación anual:

Para la información de precipitación también hay que realizaruna operación, una de las condiciones para ser zona optima es tener una precipitación anual de 1500mm a 3500mm, también puede usarse la precipitación mensual pero evaluaremos la anual para tener control de todo el año y no solo de un mes, para esto, se suman todos los rasteres de precipitación mensual para conformar la precipitación anual acumulada, la cual se obtiene de la siguiente manera:

#PRECIPITACION ANUAL
# Ruta base de los rasters mensuales
#rutaPre = "C:/Users/hilla/OneDrive/Escritorio/Semestre 5/Analisis de datos espaciales/Taller1_analisis/wc2.1_30s_prec/"

# Crear una lista con los nombres de los archivos de los 12 meses
archivos_meses = list.files(path = rutaPre, pattern = "wc2.1_2.5m_prec_.*\\.tif$", full.names = TRUE)

# Cargar todos los rasters en una lista y combinarlos
rasters_mensuales = rast(archivos_meses)

# Calcular la suma anual de precipitación
Suma_Prec = sum(rasters_mensuales)
## |---------|---------|---------|---------|=========================================                                          
# Mostrar el raster de la suma anual
plot(Suma_Prec, main = "Suma Anual de Precipitación")

PROCESAMIENTO

Clasificación de las variables:

El rrpimer paso del procesamiento consta en clasificar la temperatura y la precipitación. Luego de haber organizado la información original se toman los rangos optimos para cada variable y se clasifican como: zona apta=1 y zona no apta =0, esta clasificación aplica para ambas variables.

#Temperatura
# Clasificación de temperatura (1 si está en el rango, 0 si no)
apt_temp = (prom_Temp >= 22.5 & prom_Temp <= 28) * 1
## |---------|---------|---------|---------|=========================================                                          |---------|---------|---------|---------|=========================================                                          |---------|---------|---------|---------|=========================================                                          |---------|---------|---------|---------|=========================================                                          
# Visualización del resultado
plot(apt_temp, main = "Aptitud por Temperatura Media")

#Precipitación
# Clasificación de precipitación mensual (1 si está en el rango, 0 si no)
apt_prec = (Suma_Prec >= 1500 & Suma_Prec <= 3500) * 1
## |---------|---------|---------|---------|=========================================                                          |---------|---------|---------|---------|=========================================                                          |---------|---------|---------|---------|=========================================                                          |---------|---------|---------|---------|=========================================                                          
# Visualización del resultado
plot(apt_prec, main = "Aptitud por Precipitación Anual")

Suma de clasificaciones:

Debido a que la clasificación se divide en dos categorias: apto y no apto, se realizará una operación para hallar las zonas que son completamente aptas, o sea, que se encuentren en un rango optimo para temperatura y un rango optimo para precipitación, la operación que se utilizará es la suma, la suma de los rasteres clasificados garantiza que los que tengan el resultado mayor (2) sean zonas optimas, las zonas que obtengan como resultado 1 solo serán optimas para una variable y las zonas igual a 0 no serán una zona optima en ninguna variable.

#SUMA DE RASTERES PARA HALLAR LAS ZONAS APTAS
# Suma de ambas clasificaciones
apt_combined = apt_temp + apt_prec
## |---------|---------|---------|---------|=========================================                                          
# Visualización de la clasificación combinada
plot(apt_combined, main = "Clasificación Combinada")

Mapa de zonas optimas para cultivo de caña de azucar:

Luego de procesar la información se crea un mapa con las zonas climáticas aptas para el cultivo de caña de azucar, teniendo en cuenta el procesamiento que se venía realizando las zonas aptas son las que corresponden a la sumatoria de clasificacines igual a 2, las cuales se pueden observar en la siguiente gráfica:

# Filtrar las zonas con valor igual a 2
zonas_iguales_2 = apt_combined == 2
## |---------|---------|---------|---------|=========================================                                          
# Mostrar las zonas filtradas
plot(zonas_iguales_2, main = "Zonas con Valor 2")

Paises con zonas optimas:

Luego de tener la clasificación de zonas optimas y no optimas se crea una salida gráfica donde se observan unicamente las zonas idoneas para la siembra de caña de azucar, a continuación se observa el mapa:

# Convertir las zonas óptimas a un objeto vectorial (polígonos)
zonas_optimas_vect = as.polygons(zonas_iguales_2, values = TRUE)
## |---------|---------|---------|---------|=========================================                                          
# Filtrar solo las zonas con valor 2
zonas_optimas_vect = zonas_optimas_vect[zonas_optimas_vect$layer == 2, ]

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.

Teniendo en cuenta la información del mapa anterior, seleccionamos algunos paises que se encuentren dentro de la clasificación apta.

#INDONESIA
# Recortar el shapefile de Indonesia
Shape_Indonesia <- paises[paises$ADM0_NAME == "Indonesia", ]

# Recortar el rasterfile de Indonesia en el mapa de precipitaciones
Raster_Indonesia <- crop(zonas_iguales_2, Shape_Indonesia)

# Graficar el mapa de Indonesia con las zonas de condiciones climáticas deseadas
plot(Raster_Indonesia, legend=FALSE, main = "Zonas de Indonesia\n con las condiciones climáticas deseadas")
plot(st_geometry(Shape_Indonesia), add = TRUE, border = "blue", lwd = 2)

Como se puede observar en el mapa Indonesia es un pais que tiene una gran extensión de zona apta para el cultivo de caña de azucar, en el gráfico del país se aprecia con claridad, aun así en el mapa global se nota igualmente, la mayor parte de su territorio cuenta con las condiciones climáticas idoneas para cultivar caña de azucar, por lo tanto, es una buena opción.

#COLOMBIA
# Recortar el shapefile de Indonesia
Shape_Colombia <- paises[paises$ADM0_NAME == "Colombia", ]

# Recortar el rasterfile de Indonesia en el mapa de precipitaciones
Raster_Colombia <- crop(zonas_iguales_2, Shape_Colombia)

# Graficar el mapa de Indonesia con las zonas de condiciones climáticas deseadas
plot(Raster_Colombia, legend=FALSE, main = "Zonas de Colombia\n con las condiciones climáticas deseadas")
plot(st_geometry(Shape_Colombia), add = TRUE, border = "blue", lwd = 2)

Colombia también es un país con una gran extensión de territorio apto, no lo es tanto como Indonesia pero aun así, a simple vista es nontorio que más de la mitad desu territorio es apto para la siembra de esta planta, además de que realmente es un país con grandes producciones de cañas de azucar, por lo tanto, es un muy buen candidato para sembrarla.