Mapas de aptitud climática y similaridad para caña de azúcar
Análisis espacial global con datos climáticos de línea base
Luis Castillo
09/05/2026
1 Contexto del caso
El objetivo de esta actividad es construir una solución reproducible en RMarkdown para analizar la aptitud climática potencial de la caña de azúcar a escala global, usando información climática de línea base y herramientas de análisis geográfico.
La caña de azúcar es un cultivo tropical que requiere condiciones cálidas, suficiente disponibilidad de humedad y una estación de crecimiento prolongada. Para esta actividad se trabaja únicamente con variables climáticas, por lo cual los resultados deben interpretarse como una aptitud climática potencial.
2 Fuentes de datos y criterios técnicos
2.1 Datos climáticos
Se utiliza WorldClim versión 2.1, fuente global de
capas climáticas en formato ráster. Para mantener un equilibrio entre
detalle espacial y tiempo de procesamiento, la resolución por defecto se
define en 10 minutos de arco. Si se cuenta con un
computador de mayor capacidad, se puede cambiar a 5,
2.5 o 0.5, teniendo en cuenta que el tiempo de
descarga y procesamiento aumentará.
Variables usadas:
- Temperatura media mensual (
tavg). - Precipitación mensual (
prec).
A partir de estas capas se calculan:
- Temperatura media anual.
- Precipitación anual acumulada.
- Índices de aptitud climática por temperatura y precipitación.
- Índice combinado de aptitud climática.
- Perfiles mensuales para análisis de similaridad.
2.2 Rangos climáticos de referencia para caña de azúcar
Para la actividad se aplican rangos climáticos simplificados basados en requerimientos generales del cultivo:
| Variable | Rango marginal inferior | Rango óptimo | Rango marginal superior | Interpretación |
|---|---|---|---|---|
| Temperatura media anual | 20 a 22 °C | 22 a 30 °C | 30 a 34 °C | La temperatura limita crecimiento si es muy baja o excesivamente alta. |
| Precipitación anual | 1000 a 1500 mm | 1500 a 2500 mm | 2500 a 3500 mm | Se aproxima la disponibilidad hídrica con precipitación anual. |
la FAO reporta que la caña de azúcar alcanza crecimiento óptimo con temperaturas medias diarias entre 22 y 30 °C y que su requerimiento hídrico total puede estar entre 1500 y 2500 mm durante el ciclo de crecimiento. En este trabajo se usa la precipitación anual como aproximación climática general, aunque no reemplaza un balance hídrico completo.
2.3 Métrica de aptitud
Se utiliza una función de pertenencia gradual entre 0 y 1:
0: condición climática no apta o fuertemente restrictiva.- Entre
0y1: condición marginal o parcialmente apta. 1: condición dentro del rango óptimo.
El índice final se calcula con el criterio del factor limitante:
\[ Aptitud = \min(Aptitud_{temperatura}, Aptitud_{precipitación}) \]
Este criterio es conservador: si una zona cumple muy bien en temperatura, pero falla en precipitación, la aptitud final baja.
3 Carga y preparación de datos
# Carpeta local donde se descargarán los datos.
path_datos <- "datos_climaticos_worldclim"
dir.create(path_datos, showWarnings = FALSE, recursive = TRUE)
# Resolución WorldClim en minutos de arco.
# Opciones válidas: 10, 5, 2.5, 0.5.
# Para computadores con recursos limitados se recomienda 10.
resolucion_worldclim <- 10
# Umbral para clasificar áreas de alto potencial climático.
umbral_alto_potencial <- 0.75# Descarga de datos climáticos de línea base desde WorldClim.
# La primera ejecución puede tardar porque descarga los archivos.
tavg <- geodata::worldclim_global(
var = "tavg",
res = resolucion_worldclim,
path = path_datos
)
prec <- geodata::worldclim_global(
var = "prec",
res = resolucion_worldclim,
path = path_datos
)
# Descarga de límites administrativos globales.
paises <- geodata::world(
resolution = 3,
path = path_datos
)# Nombrar capas mensuales de forma clara.
names(tavg) <- paste0("tavg_", stringr::str_pad(1:12, width = 2, pad = "0"))
names(prec) <- paste0("prec_", stringr::str_pad(1:12, width = 2, pad = "0"))
# Algunos repositorios entregan temperatura en grados Celsius multiplicados por 10.
# Esta función detecta valores máximos inusualmente altos y ajusta automáticamente.
ajustar_temperatura <- function(r) {
maximo <- terra::global(r[[1]], fun = "max", na.rm = TRUE)[1, 1]
if (!is.na(maximo) && maximo > 60) {
r <- r / 10
}
return(r)
}
tavg <- ajustar_temperatura(tavg)
# Identificación robusta de campos en el shape global.
# geodata::world puede entregar distintos nombres de columnas según la versión
# instalada y la fuente descargada. Por eso se estandarizan nombres internos.
campos_paises <- names(paises)
candidatos_nombre <- c(
"NAME_0", "COUNTRY", "NAME", "ADMIN", "NAME_LONG", "FORMAL_EN", "SOVEREIGNT"
)
candidatos_codigo <- c(
"GID_0", "ISO_A3", "ISO3", "ADM0_A3", "SOV_A3", "WB_A3", "ISO", "GID"
)
col_nombre_pais <- candidatos_nombre[candidatos_nombre %in% campos_paises][1]
col_codigo_pais <- candidatos_codigo[candidatos_codigo %in% campos_paises][1]
if (is.na(col_nombre_pais) || length(col_nombre_pais) == 0) {
col_nombre_pais <- campos_paises[1]
}
if (is.na(col_codigo_pais) || length(col_codigo_pais) == 0) {
col_codigo_pais <- col_nombre_pais
}
# Crear campos internos estandarizados. De aquí en adelante el análisis usa
# pais_nombre, pais_codigo y pais_idx, no los nombres originales del shape.
# Se usa as.data.frame() para extraer atributos escalares del SpatVector y evitar
# que en títulos o tablas aparezcan vectores largos de países.
atributos_originales_paises <- as.data.frame(paises)
paises$pais_nombre <- as.character(atributos_originales_paises[[col_nombre_pais]])
paises$pais_codigo <- as.character(atributos_originales_paises[[col_codigo_pais]])
# Retirar Antártida si aparece en la capa global.
atributos_paises_tmp <- as.data.frame(paises)
filas_no_antartida <- !tolower(atributos_paises_tmp$pais_nombre) %in% c("antarctica", "antártida")
paises <- paises[filas_no_antartida, ]
# Índice interno después de retirar Antártida. Este índice es el que permite
# seleccionar países sin depender de códigos ISO, GID o nombres que puedan variar.
paises$pais_idx <- seq_len(nrow(paises))
# Validación informativa de los campos usados.
tibble::tibble(
campo_nombre_original = col_nombre_pais,
campo_codigo_original = col_codigo_pais,
numero_paises = nrow(paises)
) %>%
knitr::kable(caption = "Campos identificados en el shape global de países")| campo_nombre_original | campo_codigo_original | numero_paises |
|---|---|---|
| NAME_0 | GID_0 | 244 |
# Estandarización espacial preventiva.
# WorldClim trabaja en coordenadas geográficas. El shape global debe quedar en el
# mismo sistema de referencia para evitar errores en crop/mask.
crs_clima <- terra::crs(tavg, proj = TRUE)
crs_paises <- terra::crs(paises, proj = TRUE)
if (!is.na(crs_clima) && !is.na(crs_paises) && crs_clima != crs_paises) {
paises <- terra::project(paises, crs_clima)
}
# Corrección preventiva de geometrías, disponible en versiones recientes de terra.
if ("makeValid" %in% getNamespaceExports("terra")) {
paises <- terra::makeValid(paises)
}# Temperatura media anual y precipitación anual acumulada.
temp_media_anual <- mean(tavg, na.rm = TRUE)
prec_anual <- sum(prec, na.rm = TRUE)
names(temp_media_anual) <- "temperatura_media_anual_c"
names(prec_anual) <- "precipitacion_anual_mm"4 Construcción del mapa global de aptitud climática
4.1 Función de aptitud por rangos óptimos
calcular_aptitud_rango <- function(x, absoluto_min, optimo_min, optimo_max, absoluto_max) {
aptitud <- terra::ifel(
x < absoluto_min | x > absoluto_max,
0,
terra::ifel(
x >= optimo_min & x <= optimo_max,
1,
terra::ifel(
x < optimo_min,
(x - absoluto_min) / (optimo_min - absoluto_min),
(absoluto_max - x) / (absoluto_max - optimo_max)
)
)
)
aptitud <- terra::clamp(aptitud, lower = 0, upper = 1, values = TRUE)
return(aptitud)
}4.2 Cálculo de aptitud por temperatura, precipitación y aptitud combinada
apt_temp <- calcular_aptitud_rango(
temp_media_anual,
absoluto_min = 20,
optimo_min = 22,
optimo_max = 30,
absoluto_max = 34
)
apt_prec <- calcular_aptitud_rango(
prec_anual,
absoluto_min = 1000,
optimo_min = 1500,
optimo_max = 2500,
absoluto_max = 3500
)
names(apt_temp) <- "aptitud_temperatura"
names(apt_prec) <- "aptitud_precipitacion"
# Índice combinado por factor limitante.
# Nota técnica:
# Algunas versiones de terra no exportan terra::min() para combinar rasters.
# Por eso se usa terra::app(), que calcula el mínimo celda a celda de forma compatible.
aptitud_climatica <- terra::app(
c(apt_temp, apt_prec),
fun = function(x) {
if (all(is.na(x))) {
return(NA_real_)
}
min(x, na.rm = TRUE)
}
)
names(aptitud_climatica) <- "aptitud_climatica"4.3 Funciones de graficación
graficar_raster_global <- function(r, titulo, subtitulo = NULL, leyenda = "Valor") {
ggplot2::ggplot() +
tidyterra::geom_spatraster(data = r) +
tidyterra::geom_spatvector(
data = paises,
fill = NA,
color = "grey30",
linewidth = 0.12
) +
ggplot2::scale_fill_viridis_c(
name = leyenda,
limits = c(0, 1),
na.value = "transparent",
option = "C"
) +
ggplot2::coord_sf(expand = FALSE) +
ggplot2::labs(
title = titulo,
subtitle = subtitulo,
x = "Longitud",
y = "Latitud"
) +
ggplot2::theme(
legend.position = "right",
panel.grid.major = element_line(linewidth = 0.1, color = "grey85"),
plot.title = element_text(face = "bold")
)
}
graficar_variable_climatica <- function(r, titulo, subtitulo = NULL, leyenda = "Valor") {
ggplot2::ggplot() +
tidyterra::geom_spatraster(data = r) +
tidyterra::geom_spatvector(
data = paises,
fill = NA,
color = "grey30",
linewidth = 0.12
) +
ggplot2::scale_fill_viridis_c(
name = leyenda,
na.value = "transparent",
option = "D"
) +
ggplot2::coord_sf(expand = FALSE) +
ggplot2::labs(
title = titulo,
subtitle = subtitulo,
x = "Longitud",
y = "Latitud"
) +
ggplot2::theme(
legend.position = "right",
panel.grid.major = element_line(linewidth = 0.1, color = "grey85"),
plot.title = element_text(face = "bold")
)
}4.4 Mapas climáticos de entrada
graficar_variable_climatica(
temp_media_anual,
titulo = "Temperatura media anual global",
subtitulo = "Fuente climática: WorldClim 2.1, línea base 1970-2000",
leyenda = "°C"
)
este mapa permite identificar las zonas cálidas del planeta donde la temperatura media anual se aproxima a los rangos requeridos por la caña de azúcar. Sin embargo, una temperatura favorable no garantiza aptitud agrícola si la precipitación es insuficiente o excesiva.
graficar_variable_climatica(
prec_anual,
titulo = "Precipitación anual global",
subtitulo = "Suma de precipitación mensual promedio de WorldClim 2.1",
leyenda = "mm/año"
)
este mapa muestra la disponibilidad hídrica climática aproximada. Para caña de azúcar, las áreas con precipitación anual cercana a 1500-2500 mm son climáticamente favorables bajo el enfoque simplificado de esta actividad.
4.5 Mapa global de aptitud climática para caña de azúcar
graficar_raster_global(
aptitud_climatica,
titulo = "Mapa global de aptitud climática para caña de azúcar",
subtitulo = "Índice 0-1 calculado por factor limitante entre temperatura y precipitación",
leyenda = "Aptitud"
)
las zonas con valores cercanos a 1 representan áreas
donde la temperatura media anual y la precipitación anual se encuentran
dentro o cerca de los rangos óptimos definidos para el cultivo. Las
zonas con valores bajos no necesariamente son imposibles para la
producción, pero requerirían condiciones no evaluadas aquí, como riego,
manejo agronómico intensivo, variedades específicas o infraestructura
adicional.
graficar_raster_global(
apt_temp,
titulo = "Componente de aptitud por temperatura",
subtitulo = "Rango óptimo usado: 22 a 30 °C",
leyenda = "Aptitud"
)
graficar_raster_global(
apt_prec,
titulo = "Componente de aptitud por precipitación",
subtitulo = "Rango óptimo usado: 1500 a 2500 mm/año",
leyenda = "Aptitud"
)
comparar estos dos mapas permite reconocer cuál variable actúa como limitante principal. Una zona cálida puede tener baja aptitud si es muy seca; de igual manera, una zona húmeda puede ser limitada si su temperatura media anual es baja.
5 Identificación de países con alto potencial climático
5.1 Estimación de área de alto potencial por país
# Área de cada celda en km2.
area_celdas_km2 <- terra::cellSize(aptitud_climatica, unit = "km")
# Solo se suma el área donde la aptitud supera el umbral definido.
area_alta <- terra::ifel(aptitud_climatica >= umbral_alto_potencial, area_celdas_km2, 0)
names(area_alta) <- "area_alta_km2"
area_extraida <- terra::extract(
area_alta,
paises,
fun = sum,
na.rm = TRUE
)
# La salida de terra::extract incluye un campo ID que corresponde a la posición
# del país dentro del objeto SpatVector. Se usa ese ID para evitar errores de
# coincidencia por nombres o códigos de país.
area_extraida_df <- as.data.frame(area_extraida)
columna_area <- setdiff(names(area_extraida_df), "ID")[1]
atributos_paises_df <- as.data.frame(paises)
atributos_paises <- tibble::tibble(
pais_idx = as.integer(atributos_paises_df$pais_idx),
pais_codigo = as.character(atributos_paises_df$pais_codigo),
pais_nombre = as.character(atributos_paises_df$pais_nombre)
)
area_paises <- tibble::tibble(
pais_idx = as.integer(area_extraida_df$ID),
area_alta_km2 = as.numeric(area_extraida_df[[columna_area]])
) %>%
dplyr::left_join(atributos_paises, by = "pais_idx") %>%
dplyr::mutate(
area_alta_km2 = dplyr::if_else(is.na(area_alta_km2), 0, area_alta_km2)
) %>%
dplyr::arrange(dplyr::desc(area_alta_km2))
# Selección principal: países con área de alto potencial.
top_paises <- area_paises %>%
dplyr::filter(area_alta_km2 > 0) %>%
dplyr::slice_head(n = 3) %>%
dplyr::select(pais_idx, pais_codigo, pais_nombre, area_alta_km2)
# Fallback técnico: si el umbral elegido no genera áreas > 0, se toman los tres
# países con mayor valor calculado para que el documento siga ejecutándose y se
# deje trazabilidad de la situación en la interpretación.
if (nrow(top_paises) == 0) {
top_paises <- area_paises %>%
dplyr::slice_head(n = 3) %>%
dplyr::select(pais_idx, pais_codigo, pais_nombre, area_alta_km2)
mensaje_top_paises <- paste0(
"Con el umbral de alto potencial definido (", umbral_alto_potencial,
") no se encontraron países con área > 0. Se muestran los tres países con ",
"mayor valor relativo calculado para continuar el análisis."
)
} else {
mensaje_top_paises <- paste0(
"Países con mayor área estimada de alto potencial climático. Umbral de aptitud >= ",
umbral_alto_potencial
)
}
# Se construye una versión de presentación para que la tabla sea clara y
# no muestre campos técnicos usados internamente para el recorte espacial.
total_area_alta_global <- sum(area_paises$area_alta_km2, na.rm = TRUE)
top_paises_tabla <- top_paises %>%
dplyr::mutate(
puesto = dplyr::row_number(),
participacion_pct = if (is.finite(total_area_alta_global) && total_area_alta_global > 0) {
(area_alta_km2 / total_area_alta_global) * 100
} else {
rep(0, dplyr::n())
},
area_formato = format(
round(area_alta_km2, 2),
big.mark = ".",
decimal.mark = ",",
nsmall = 2,
scientific = FALSE
),
participacion_formato = paste0(
format(
round(participacion_pct, 2),
decimal.mark = ",",
nsmall = 2,
scientific = FALSE
),
"%"
)
) %>%
dplyr::transmute(
`Puesto` = puesto,
`País` = pais_nombre,
`Código` = pais_codigo,
`Área de alto potencial climático (km²)` = area_formato,
`Participación sobre el área global de alto potencial` = participacion_formato
)
knitr::kable(
top_paises_tabla,
format = "html",
escape = FALSE,
align = c("c", "l", "c", "r", "r"),
caption = mensaje_top_paises
) %>%
kableExtra::kable_styling(
bootstrap_options = c("striped", "hover", "condensed", "responsive"),
full_width = FALSE,
position = "center",
font_size = 13
) %>%
kableExtra::row_spec(0, bold = TRUE) %>%
kableExtra::column_spec(1, width = "8em") %>%
kableExtra::column_spec(2, width = "16em") %>%
kableExtra::column_spec(3, width = "8em") %>%
kableExtra::column_spec(4, width = "15em") %>%
kableExtra::column_spec(5, width = "18em")| Puesto | País | Código | Área de alto potencial climático (km²) | Participación sobre el área global de alto potencial |
|---|---|---|---|---|
| 1 | Brazil | BRA | 5.243.249,95 | 33,65% |
| 2 | Democratic Republic of the Congo | COD | 1.724.172,60 | 11,07% |
| 3 | Indonesia | IDN | 876.586,13 | 5,63% |
la tabla anterior identifica los países que concentran mayor superficie climáticamente favorable bajo los criterios definidos
5.2 Recorte espacial para países con alto potencial
# Selección robusta de países top.
# Se seleccionan las filas del SpatVector usando el campo interno pais_idx.
# Esto evita que los títulos se construyan con listas largas de países o que se
# recorten geometrías no esperadas por una lectura ambigua de atributos.
indices_top <- as.integer(top_paises$pais_idx)
indices_top <- indices_top[!is.na(indices_top)]
atributos_paises_recorte <- as.data.frame(paises)
filas_top <- which(as.integer(atributos_paises_recorte$pais_idx) %in% indices_top)
if (length(filas_top) == 0) {
stop("No se encontraron países para recortar. Revise la tabla top_paises y el índice interno pais_idx.")
}
paises_top <- paises[filas_top, ]
nombres_top <- atributos_paises_recorte$pais_nombre[filas_top]
nombres_top_texto <- stringr::str_wrap(paste(nombres_top, collapse = ", "), width = 95)
# Función auxiliar para recortar un raster con un vector de forma segura.
# Primero recorta por la extensión rectangular del vector y luego aplica máscara
# con la geometría real. Esto evita el error '[crop] invalid extent' cuando el
# objeto vectorial llega con geometrías problemáticas o CRS no equivalente.
recortar_raster_por_vector <- function(raster_obj, vector_obj, nombre_vector = "vector") {
if (!inherits(vector_obj, "SpatVector")) {
stop(nombre_vector, " debe ser un objeto SpatVector de terra.")
}
if (nrow(vector_obj) == 0) {
stop(nombre_vector, " no contiene geometrías para recortar.")
}
crs_raster <- terra::crs(raster_obj, proj = TRUE)
crs_vector <- terra::crs(vector_obj, proj = TRUE)
if (!is.na(crs_raster) && !is.na(crs_vector) && crs_raster != crs_vector) {
vector_obj <- terra::project(vector_obj, crs_raster)
}
if ("makeValid" %in% getNamespaceExports("terra")) {
vector_obj <- terra::makeValid(vector_obj)
}
extension_vector <- terra::ext(vector_obj)
valores_extension <- as.vector(extension_vector)
if (
length(valores_extension) != 4 ||
any(!is.finite(valores_extension)) ||
valores_extension[1] >= valores_extension[2] ||
valores_extension[3] >= valores_extension[4]
) {
stop("La extensión espacial de ", nombre_vector, " no es válida para realizar el recorte.")
}
raster_recortado <- terra::crop(raster_obj, extension_vector, snap = "out")
raster_mascara <- terra::mask(raster_recortado, vector_obj)
return(list(raster = raster_mascara, vector = vector_obj))
}
recorte_top <- recortar_raster_por_vector(
aptitud_climatica,
paises_top,
nombre_vector = "paises_top"
)
aptitud_top <- recorte_top$raster
paises_top <- recorte_top$vector
mapa_top_paises <- ggplot2::ggplot() +
tidyterra::geom_spatraster(data = aptitud_top) +
tidyterra::geom_spatvector(
data = paises_top,
fill = NA,
color = "grey15",
linewidth = 0.35
) +
ggplot2::scale_fill_viridis_c(
name = "Aptitud",
limits = c(0, 1),
na.value = "transparent",
option = "C"
) +
ggplot2::coord_sf(expand = FALSE) +
ggplot2::labs(
title = "Recorte de países con mayor área de alto potencial climático",
subtitle = paste0("Países seleccionados:\n", nombres_top_texto),
x = "Longitud",
y = "Latitud"
) +
ggplot2::theme(
legend.position = "right",
plot.title = element_text(face = "bold")
)
mapa_top_paises
el recorte permite observar con mayor detalle la distribución interna de la aptitud en los países seleccionados. Las zonas de mayor aptitud suelen concentrarse en regiones tropicales y subtropicales con condiciones cálidas y disponibilidad hídrica suficiente.
graficar_pais_individual <- function(fila_pais_i) {
fila_pais_i <- as.integer(fila_pais_i)
pais_i <- paises[fila_pais_i, ]
atributos_i <- as.data.frame(pais_i)
nombre_i <- as.character(atributos_i$pais_nombre[1])
codigo_i <- as.character(atributos_i$pais_codigo[1])
if (is.na(nombre_i) || nombre_i == "") {
nombre_i <- paste("País", fila_pais_i)
}
subtitulo_i <- if (!is.na(codigo_i) && codigo_i != "" && codigo_i != nombre_i) {
paste0("Recorte nacional del índice climático 0-1 | Código: ", codigo_i)
} else {
"Recorte nacional del índice climático 0-1"
}
recorte_i <- recortar_raster_por_vector(
aptitud_climatica,
pais_i,
nombre_vector = paste("país", nombre_i)
)
r_i <- recorte_i$raster
pais_i <- recorte_i$vector
g <- ggplot2::ggplot() +
tidyterra::geom_spatraster(data = r_i) +
tidyterra::geom_spatvector(
data = pais_i,
fill = NA,
color = "grey15",
linewidth = 0.35
) +
ggplot2::scale_fill_viridis_c(
name = "Aptitud",
limits = c(0, 1),
na.value = "transparent",
option = "C"
) +
ggplot2::coord_sf(expand = FALSE) +
ggplot2::labs(
title = stringr::str_wrap(paste("Aptitud climática para caña de azúcar en", nombre_i), width = 70),
subtitle = subtitulo_i,
x = "Longitud",
y = "Latitud"
) +
ggplot2::theme(
legend.position = "right",
plot.title = element_text(face = "bold", size = 13),
plot.subtitle = element_text(size = 10),
plot.margin = ggplot2::margin(t = 10, r = 10, b = 10, l = 10)
)
print(g)
}
purrr::walk(filas_top, graficar_pais_individual)
6 Puntos de análisis en el Valle del Cauca
Para esta sección se seleccionan tres puntos de referencia dentro del Valle del Cauca. Las coordenadas son aproximadas y pueden verificarse en Google Maps. Se escogieron ubicaciones asociadas a zonas representativas del valle geográfico del río Cauca, donde históricamente existe actividad agrícola relacionada con caña de azúcar.
puntos_valle <- tibble::tibble(
sitio = c("Palmira", "Buga", "Tuluá"),
longitud = c(-76.3036, -76.3020, -76.1954),
latitud = c(3.5394, 3.9009, 4.0847)
)
puntos_vect <- terra::vect(
puntos_valle,
geom = c("longitud", "latitud"),
crs = "EPSG:4326"
)
knitr::kable(
puntos_valle,
digits = 4,
caption = "Puntos seleccionados en el Valle del Cauca para extracción climática"
)| sitio | longitud | latitud |
|---|---|---|
| Palmira | -76.3036 | 3.5394 |
| Buga | -76.3020 | 3.9009 |
| Tuluá | -76.1954 | 4.0847 |
# Recorte aproximado para Valle del Cauca y alrededores.
ext_valle <- terra::ext(-77.2, -75.7, 3.0, 4.7)
apt_valle <- terra::crop(aptitud_climatica, ext_valle)
paises_valle <- terra::crop(paises, ext_valle)
mapa_puntos <- ggplot2::ggplot() +
tidyterra::geom_spatraster(data = apt_valle) +
tidyterra::geom_spatvector(
data = paises_valle,
fill = NA,
color = "grey30",
linewidth = 0.25
) +
tidyterra::geom_spatvector(
data = puntos_vect,
color = "black",
fill = "white",
shape = 21,
size = 3,
stroke = 1
) +
ggplot2::scale_fill_viridis_c(
name = "Aptitud",
limits = c(0, 1),
na.value = "transparent",
option = "C"
) +
ggplot2::coord_sf(expand = FALSE) +
ggplot2::labs(
title = "Puntos climáticos seleccionados en el Valle del Cauca",
subtitle = "Base cartográfica sobre el índice de aptitud climática",
x = "Longitud",
y = "Latitud"
) +
ggplot2::theme(
legend.position = "right",
plot.title = element_text(face = "bold")
)
mapa_puntos
el mapa ubica los puntos seleccionados en el contexto regional. Estos puntos se usarán como perfiles climáticos de referencia para comparar otras regiones del mundo por similaridad climática.
7 Extracción de información climática y series de tiempo
7.1 Extracción mensual de temperatura y precipitación
extraer_series <- function(r, puntos, variable) {
datos <- terra::extract(r, puntos) %>%
tibble::as_tibble() %>%
dplyr::select(-ID) %>%
dplyr::bind_cols(puntos_valle %>% dplyr::select(sitio), .) %>%
tidyr::pivot_longer(
cols = -sitio,
names_to = "capa",
values_to = "valor"
) %>%
dplyr::mutate(
mes = as.integer(stringr::str_extract(capa, "[0-9]+$")),
variable = variable
) %>%
dplyr::select(sitio, variable, mes, valor)
return(datos)
}
serie_temp <- extraer_series(tavg, puntos_vect, "Temperatura media") %>%
dplyr::rename(temperatura_c = valor)
serie_prec <- extraer_series(prec, puntos_vect, "Precipitación") %>%
dplyr::rename(precipitacion_mm = valor)
resumen_puntos <- serie_temp %>%
dplyr::group_by(sitio) %>%
dplyr::summarise(
temperatura_media_anual_c = mean(temperatura_c, na.rm = TRUE),
.groups = "drop"
) %>%
dplyr::left_join(
serie_prec %>%
dplyr::group_by(sitio) %>%
dplyr::summarise(
precipitacion_anual_mm = sum(precipitacion_mm, na.rm = TRUE),
.groups = "drop"
),
by = "sitio"
)
knitr::kable(
resumen_puntos,
digits = 2,
caption = "Resumen climático anual de los puntos seleccionados"
)| sitio | temperatura_media_anual_c | precipitacion_anual_mm |
|---|---|---|
| Buga | 21.72 | 1600 |
| Palmira | 22.41 | 1484 |
| Tuluá | 23.59 | 1551 |
7.2 Serie mensual de temperatura
grafico_temp <- ggplot2::ggplot(
serie_temp,
ggplot2::aes(x = mes, y = temperatura_c, group = sitio, color = sitio)
) +
ggplot2::geom_line(linewidth = 1) +
ggplot2::geom_point(size = 2) +
ggplot2::scale_x_continuous(breaks = 1:12) +
ggplot2::labs(
title = "Serie mensual de temperatura media en puntos del Valle del Cauca",
subtitle = "Valores promedio mensuales de WorldClim 2.1",
x = "Mes",
y = "Temperatura media (°C)",
color = "Sitio"
) +
ggplot2::theme(
legend.position = "bottom",
plot.title = element_text(face = "bold")
)
grafico_temp
la serie de temperatura permite evaluar la estabilidad térmica mensual en los puntos seleccionados. En zonas tropicales como el Valle del Cauca se espera una variación mensual relativamente baja frente a regiones templadas, lo cual puede favorecer cultivos de ciclo largo como la caña de azúcar.
7.3 Serie mensual de precipitación
grafico_prec <- ggplot2::ggplot(
serie_prec,
ggplot2::aes(x = mes, y = precipitacion_mm, group = sitio, color = sitio)
) +
ggplot2::geom_line(linewidth = 1) +
ggplot2::geom_point(size = 2) +
ggplot2::scale_x_continuous(breaks = 1:12) +
ggplot2::labs(
title = "Serie mensual de precipitación en puntos del Valle del Cauca",
subtitle = "Valores promedio mensuales de WorldClim 2.1",
x = "Mes",
y = "Precipitación mensual (mm)",
color = "Sitio"
) +
ggplot2::theme(
legend.position = "bottom",
plot.title = element_text(face = "bold")
)
grafico_prec
la serie de precipitación permite identificar meses más húmedos o secos. Para la caña de azúcar, una distribución adecuada de humedad durante el periodo vegetativo es importante; sin embargo, para maduración y cosecha suelen ser deseables condiciones relativamente más secas.
8 Mapas globales de similaridad climática
8.1 Enfoque metodológico
Para medir similaridad climática se compara el perfil mensual de cada punto del Valle del Cauca contra el perfil mensual de cada celda del mapa global. Se usan 24 variables:
- 12 valores de temperatura media mensual.
- 12 valores de precipitación mensual.
Como temperatura y precipitación tienen unidades y escalas distintas, primero se estandarizan las variables globalmente mediante puntajes z:
\[ Z = \frac{x - \mu}{\sigma} \]
Luego se calcula la distancia euclidiana entre cada celda global y cada punto de referencia:
\[ D = \sqrt{\sum_{i=1}^{24}(x_i - r_i)^2} \]
Finalmente, la distancia se transforma a un índice de similaridad entre 0 y 1:
\[ Similaridad = \frac{1}{1 + D} \]
Valores cercanos a 1 indican mayor parecido climático con el punto de referencia.
8.2 Estandarización de variables climáticas
clima_24 <- c(tavg, prec)
medias_globales <- terra::global(clima_24, fun = "mean", na.rm = TRUE)[, 1]
desv_globales <- terra::global(clima_24, fun = "sd", na.rm = TRUE)[, 1]
clima_24_std <- clima_24
for (i in seq_len(terra::nlyr(clima_24))) {
clima_24_std[[i]] <- (clima_24[[i]] - medias_globales[i]) / desv_globales[i]
}
names(clima_24_std) <- names(clima_24)8.3 Cálculo de distancia euclidiana y similaridad
valores_puntos_std <- terra::extract(clima_24_std, puntos_vect) %>%
tibble::as_tibble() %>%
dplyr::select(-ID)
calcular_distancia_euclidiana <- function(stack_std, valores_referencia) {
capas_cuadradas <- vector("list", terra::nlyr(stack_std))
for (i in seq_len(terra::nlyr(stack_std))) {
capas_cuadradas[[i]] <- (stack_std[[i]] - valores_referencia[i])^2
}
suma_cuadrados <- sum(do.call(c, capas_cuadradas))
distancia <- sqrt(suma_cuadrados)
return(distancia)
}
similaridad_lista <- list()
for (i in seq_len(nrow(puntos_valle))) {
distancia_i <- calcular_distancia_euclidiana(
clima_24_std,
as.numeric(valores_puntos_std[i, ])
)
similaridad_i <- 1 / (1 + distancia_i)
names(similaridad_i) <- paste0("similaridad_", puntos_valle$sitio[i])
similaridad_lista[[i]] <- similaridad_i
}
similaridad_stack <- do.call(c, similaridad_lista)8.4 Mapas de similaridad global
graficar_similaridad <- function(r, sitio) {
ggplot2::ggplot() +
tidyterra::geom_spatraster(data = r) +
tidyterra::geom_spatvector(
data = paises,
fill = NA,
color = "grey30",
linewidth = 0.12
) +
ggplot2::scale_fill_viridis_c(
name = "Similaridad",
limits = c(0, 1),
na.value = "transparent",
option = "B"
) +
ggplot2::coord_sf(expand = FALSE) +
ggplot2::labs(
title = paste("Similaridad climática global respecto a", sitio),
subtitle = "Distancia euclidiana sobre temperatura y precipitación mensual estandarizadas",
x = "Longitud",
y = "Latitud"
) +
ggplot2::theme(
legend.position = "right",
plot.title = element_text(face = "bold")
)
}
for (i in seq_len(terra::nlyr(similaridad_stack))) {
print(graficar_similaridad(similaridad_stack[[i]], puntos_valle$sitio[i]))
}
los mapas de similaridad no identifican directamente zonas óptimas para caña de azúcar; identifican zonas cuyo comportamiento climático mensual se parece al de los puntos seleccionados del Valle del Cauca. Por tanto, pueden resaltar regiones climáticamente análogas incluso si no cumplen exactamente todos los rangos óptimos definidos en la primera aproximación.
9 Comparación entre aptitud por rangos y similaridad climática
# Similaridad máxima entre los tres puntos seleccionados.
# Se calcula con terra::app() para evitar incompatibilidades de métodos
# entre versiones de terra al usar max() directamente sobre un SpatRaster.
similaridad_maxima <- terra::app(
similaridad_stack,
fun = function(x) {
if (all(is.na(x))) {
return(NA_real_)
}
max(x, na.rm = TRUE)
}
)
names(similaridad_maxima) <- "similaridad_maxima_valle"
# Reescalamiento visual para comparación conjunta.
comparacion_stack <- c(aptitud_climatica, similaridad_maxima)
names(comparacion_stack) <- c("Aptitud por rangos", "Similaridad con Valle del Cauca")
ggplot2::ggplot() +
tidyterra::geom_spatraster(data = comparacion_stack) +
tidyterra::geom_spatvector(
data = paises,
fill = NA,
color = "grey30",
linewidth = 0.10
) +
ggplot2::facet_wrap(~lyr, ncol = 1) +
ggplot2::scale_fill_viridis_c(
name = "Índice",
limits = c(0, 1),
na.value = "transparent",
option = "C"
) +
ggplot2::coord_sf(expand = FALSE) +
ggplot2::labs(
title = "Comparación entre aptitud climática y similaridad climática",
subtitle = "Aptitud por rangos óptimos vs. similaridad mensual con puntos del Valle del Cauca",
x = "Longitud",
y = "Latitud"
) +
ggplot2::theme(
legend.position = "right",
plot.title = element_text(face = "bold"),
strip.text = element_text(face = "bold")
)
9.1 Lectura comparativa
La aproximación de aptitud por rangos óptimos responde a la pregunta: ¿qué zonas cumplen condiciones climáticas generales apropiadas para la caña de azúcar según temperatura y precipitación anual?
La aproximación de similaridad climática responde a una pregunta diferente: ¿qué zonas del mundo se parecen climáticamente a los puntos seleccionados del Valle del Cauca considerando la dinámica mensual de temperatura y precipitación?
Por lo tanto:
- Una zona puede tener alta aptitud por rangos y baja similaridad con el Valle del Cauca si cumple los rangos anuales, pero tiene una distribución mensual distinta.
- Una zona puede tener alta similaridad con el Valle del Cauca, pero no necesariamente ser la mejor zona productiva, porque la similaridad no evalúa suelos, pendiente, manejo agronómico ni infraestructura.
- La aptitud por rangos es útil para una evaluación agroclimática general.
- La similaridad es útil para buscar zonas análogas a un territorio conocido.
10 Conclusiones
La metodología permitió construir un índice global de aptitud climática para caña de azúcar usando variables de temperatura media anual y precipitación anual derivadas de WorldClim.
El uso de funciones de aptitud gradual evita clasificaciones rígidas y permite representar zonas marginales, intermedias y óptimas.
El índice combinado por factor limitante es adecuado para este caso porque reconoce que temperatura y precipitación son condiciones simultáneamente necesarias para el desarrollo del cultivo.
La identificación de países con mayor área de alto potencial permite pasar de una lectura global a una lectura territorial más específica, útil para análisis comparativos entre regiones.
La extracción de series climáticas en puntos del Valle del Cauca permite construir perfiles mensuales de referencia, lo cual enriquece el análisis porque no depende únicamente de promedios anuales.
Los mapas de similaridad climática permiten identificar regiones del mundo con comportamiento mensual parecido al Valle del Cauca. Esta aproximación es complementaria a la aptitud por rangos y puede ser útil para estudios de analogía climática, transferencia de conocimiento agronómico o búsqueda de zonas comparables.
La comparación entre ambas aproximaciones evidencia que la aptitud climática y la similaridad climática no son equivalentes. La primera evalúa cumplimiento de requerimientos del cultivo; la segunda evalúa parecido climático con un territorio de referencia.
Como limitación principal, este análisis no incorpora variables edáficas, topográficas, disponibilidad real de riego, restricciones ambientales, cobertura del suelo, productividad histórica ni condiciones socioeconómicas. Por tanto, los mapas deben interpretarse como una primera aproximación climática y no como una zonificación agrícola definitiva.
11 Referencias
- FAO. Crop information: Sugarcane. Food and Agriculture Organization of the United Nations.
- WorldClim. Bioclimatic variables and historical climate data, versión 2.1.
- Hijmans, R. J., Cameron, S. E., Parra, J. L., Jones, P. G., & Jarvis, A. (2005). Very high resolution interpolated climate surfaces for global land areas. International Journal of Climatology, 25(15), 1965-1978.
- Paquetes de R usados:
terra,geodata,ggplot2,tidyterra,dplyr,tidyr,purrr.