# ============================================================
# NORMALIZACION GENERAL DE TABLAS PARA HTML
# ============================================================

normalizar_tabla_utf8 <- function(tabla) {

  tabla <- as.data.frame(
    tabla,
    stringsAsFactors = FALSE,
    check.names = FALSE
  )

  # Normalizar nombres de columnas
  nombres_originales <- names(tabla)

  nombres_utf8 <- iconv(
    nombres_originales,
    from = "",
    to = "UTF-8",
    sub = "byte"
  )

  # Si el sistema no reconoce la codificación de origen,
  # intentar Windows-1252 únicamente sobre los nombres fallidos.
  nombres_fallidos <- is.na(nombres_utf8)

  if (any(nombres_fallidos)) {
    nombres_utf8[nombres_fallidos] <- iconv(
      nombres_originales[nombres_fallidos],
      from = "windows-1252",
      to = "UTF-8",
      sub = "byte"
    )
  }

  nombres_utf8 <- make.unique(
    nombres_utf8,
    sep = "_"
  )

  names(tabla) <- nombres_utf8

  # Normalizar columnas de texto
  columnas_texto <- vapply(
    tabla,
    function(x) {
      is.character(x) || is.factor(x)
    },
    logical(1)
  )

  tabla[columnas_texto] <- lapply(
    tabla[columnas_texto],
    function(x) {

      x <- as.character(x)

      convertido <- iconv(
        x,
        from = "",
        to = "UTF-8",
        sub = "byte"
      )

      fallidos <- is.na(convertido) & !is.na(x)

      if (any(fallidos)) {
        convertido[fallidos] <- iconv(
          x[fallidos],
          from = "windows-1252",
          to = "UTF-8",
          sub = "byte"
        )
      }

      convertido
    }
  )

  tabla
}

Informe 2 RMarkdown: aptitud climática global para caña de azúcar

Universidad Pontificia Javeriana

Presentado por: Dainer Orlando Cruz Chate

Trabajo presentado al profesor: David Arango Londoño

Curso: Análisis de Datos Geográficos y Espaciales

Unidad: 2

Fecha: 27/10/2026

1 Planteamiento del caso

Un agricultor desea identificar zonas del mundo con condiciones climáticas favorables para la caña de azúcar. Se aplican dos aproximaciones complementarias:

1. Aptitud por umbrales: temperatura media anual entre 22,5 y 28 °C y precipitación anual entre 1.500 y 3.500 mm.
2. Analogía climática: similitud entre el ciclo mensual de temperatura y precipitación de tres sitios del Valle del Cauca y las celdas climáticas del mundo.

Alcance: el resultado representa aptitud climática de línea base, no una predicción directa de rendimiento en kg/ha. Para estimar rendimiento potencial se requieren, además, radiación, suelo, relieve, variedad, duración del ciclo, manejo, riego y un modelo de cultivo calibrado, por ejemplo DSSAT-CANEGRO, APSIM-Sugar o AquaCrop.

2 Fuentes, criterios y método

Se emplea WorldClim 2.1, climatología de línea base 1970–2000. WorldClim ofrece capas mensuales de temperatura media y precipitación a varias resoluciones. Para que el informe pueda ejecutarse en un computador convencional se utiliza resolución de 10 minutos de arco; puede cambiarse a 5 o 2,5 minutos, con mayor costo de memoria y tiempo.

Los rangos de aptitud usados en este ejercicio proceden de las instrucciones del caso. Como referencia agronómica adicional, FAO señala que la producción de caña depende de clima, suelo, variedad y manejo, y reporta condiciones cálidas y pluviometría elevada como favorables.

Fuentes principales

• WorldClim 2.1: https://www.worldclim.org/data/worldclim21.html
• Fick, S. E. y Hijmans, R. J. (2017), WorldClim 2: https://doi.org/10.1002/joc.5086
• FAO, información general sobre caña de azúcar: https://www.fao.org/4/x5648e/x5648e0e.htm

2.1 Herramientas espaciales y compatibilidad

El flujo de trabajo integra raster para compatibilidad con las estructuras ráster solicitadas, leaflet para mapas interactivos, tmaptools para utilidades cartográficas y geocodificación, y sf/terra para las operaciones espaciales actuales. El paquete rgdal se comprueba en el entorno por exigencia del caso estudiado, pero no se fuerza su instalación porque fue retirado de CRAN; las transformaciones y lecturas se realizan con sf y terra para asegurar compatibilidad con versiones recientes de R.

# paquetes_caso <- data.frame(
#  Paquete = c("tmaptools", "leaflet", "raster", "rgdal"),
#  Función_en_el_informe = c("Utilidades cartográficas y geocodificación",
#                            "Visualización cartográfica interactiva",
#                            "Compatibilidad y manejo de capas ráster",
#                            "Compatibilidad heredada; sustituido operativamente por sf/terra"),
#  Disponible = c(requireNamespace("tmaptools", quietly = TRUE),
#                 requireNamespace("leaflet", quietly = TRUE),
#                 requireNamespace("raster", quietly = TRUE),
#                 rgdal_disponible)
# )
# knitr::kable(paquetes_caso, caption = "Paquetes espaciales contemplados") |>
#  kableExtra::kable_styling(bootstrap_options = c("striped", "hover"), full_width = TRUE)

# Tabla de paquetes espaciales utilizados o contemplados
paquetes_caso <- data.frame(
  Paquete = c(
    "tmaptools",
    "leaflet",
    "raster",
    "rgdal"
  ),
  Funcion_en_el_informe = c(
    "Utilidades cartograficas y geocodificacion",
    "Visualizacion cartografica interactiva",
    "Compatibilidad y manejo de capas raster",
    "Compatibilidad heredada; sustituido operativamente por sf/terra"
  ),
  Disponible = c(
    requireNamespace("tmaptools", quietly = TRUE),
    requireNamespace("leaflet", quietly = TRUE),
    requireNamespace("raster", quietly = TRUE),
    isTRUE(rgdal_disponible)
  ),
  stringsAsFactors = FALSE,
  check.names = FALSE
)

# Etiquetas comprensibles para la disponibilidad
paquetes_caso$Disponible <- ifelse(
  paquetes_caso$Disponible,
  "Si",
  "No"
)

# Convertir explícitamente el contenido a UTF-8
paquetes_caso[] <- lapply(
  paquetes_caso,
  function(x) {
    if (is.character(x)) {
      iconv(
        x,
        from = "",
        to = "UTF-8",
        sub = "byte"
      )
    } else {
      x
    }
  }
)

paquetes_caso <- normalizar_tabla_utf8(paquetes_caso)

# Generar la tabla HTML

tabla_paquetes_caso <- knitr::kable(
  paquetes_caso,
  format = "html",
  col.names = c(
    "Paquete",
    "Funcion en el informe",
    "Disponible"
  ),
  caption = "Paquetes espaciales contemplados",
  align = c("l", "l", "c"),
  escape = TRUE
)

kableExtra::kable_styling(
  tabla_paquetes_caso,
  bootstrap_options = c(
    "striped",
    "hover",
    "condensed",
    "responsive"
  ),
  full_width = TRUE,
  position = "center"
)

Paquetes espaciales contemplados

Paquete	Funcion en el informe	Disponible
tmaptools	Utilidades cartograficas y geocodificacion	Si
leaflet	Visualizacion cartografica interactiva	Si
raster	Compatibilidad y manejo de capas raster	Si
rgdal	Compatibilidad heredada; sustituido operativamente por sf/terra	No

2.2 Descarga y preparación de datos

# Climatologías mensuales globales. geodata almacena los archivos en ruta_datos.
tmed_m <- geodata::worldclim_global(var = "tavg", res = 10, path = ruta_datos)
prec_m <- geodata::worldclim_global(var = "prec", res = 10, path = ruta_datos)

names(tmed_m) <- paste0("tmed_", sprintf("%02d", 1:12))
names(prec_m) <- paste0("prec_", sprintf("%02d", 1:12))

# Verificación automática de unidades de temperatura.
# En algunas distribuciones raster la temperatura puede estar escalada por 10.
media_prueba <- global(tmed_m[[1]], "mean", na.rm = TRUE)[1, 1]
if (is.finite(media_prueba) && media_prueba > 60) tmed_m <- tmed_m / 10

# Agregados anuales
tmed_anual <- mean(tmed_m, na.rm = TRUE)
prec_anual <- sum(prec_m, na.rm = TRUE)
names(tmed_anual) <- "temperatura_media_C"
names(prec_anual) <- "precipitacion_anual_mm"

# Conversión explícita al formato RasterLayer solicitado en la rúbrica.
# Se conserva terra para cálculos robustos y raster para interoperabilidad.
tmed_anual_raster <- raster::raster(tmed_anual)
prec_anual_raster <- raster::raster(prec_anual)

# Límites nacionales, provistos por geodata/GADM
paises <- geodata::world(path = ruta_datos)

# Control básico
stopifnot(nlyr(tmed_m) == 12, nlyr(prec_m) == 12)
resumen_datos <- data.frame(
  Característica = c("Objeto", "Filas", "Columnas", "Capas climáticas",
                     "Resolución espacial", "Extensión", "Referencia espacial",
                     "Periodo climático"),
  Valor = c("Ráster multibanda global", terra::nrow(tmed_m), terra::ncol(tmed_m),
            terra::nlyr(tmed_m) + terra::nlyr(prec_m),
            paste(round(terra::res(tmed_m), 6), collapse = " × "),
            paste(round(as.vector(terra::ext(tmed_m)), 2), collapse = ", "),
            "WGS 84 (EPSG:4326)", "1970–2000")
)

resumen_datos <- normalizar_tabla_utf8(resumen_datos)
knitr::kable(resumen_datos, caption = "Características de los datos climáticos utilizados",
             align = c("l", "l")) |>
  kableExtra::kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                            full_width = TRUE)

Características de los datos climáticos utilizados

Característica	Valor
Objeto	Ráster multibanda global
Filas	1080
Columnas	2160
Capas climáticas	24
Resolución espacial	0.166667 × 0.166667
Extensión	-180, 180, -90, 90
Referencia espacial	WGS 84 (EPSG:4326)
Periodo climático	1970–2000

3 Aptitud climática por rangos óptimos

Se calcula primero la adecuación de cada variable en escala continua de 0 a 1 mediante una función trapezoidal con tolerancia exterior. Después se obtiene un índice conservador como el mínimo entre temperatura y precipitación. Así, una limitación fuerte en una variable no queda compensada por la otra.

# VERSIÓN 5: cálculo mediante álgebra raster nativa, sin funciones personalizadas.
# Los objetos espaciales se reconstruyen en cada renderizado y no se almacenan en caché.

# Adecuación térmica trapezoidal, calculada directamente con álgebra raster.
apt_temp <- terra::ifel(
  tmed_anual < 18 | tmed_anual > 35, 0,
  terra::ifel(
    tmed_anual < 22.5,
    (tmed_anual - 18) / (22.5 - 18),
    terra::ifel(
      tmed_anual <= 28,
      1,
      (35 - tmed_anual) / (35 - 28)
    )
  )
)

# Adecuación de precipitación trapezoidal.
apt_prec <- terra::ifel(
  prec_anual < 750 | prec_anual > 5000, 0,
  terra::ifel(
    prec_anual < 1500,
    (prec_anual - 750) / (1500 - 750),
    terra::ifel(
      prec_anual <= 3500,
      1,
      (5000 - prec_anual) / (5000 - 3500)
    )
  )
)

# El operador conservador toma la menor adecuación de las dos variables.
apt_climatica <- terra::ifel(apt_temp <= apt_prec, apt_temp, apt_prec)
names(apt_temp) <- "aptitud_temperatura"
names(apt_prec) <- "aptitud_precipitacion"
names(apt_climatica) <- "aptitud_climatica"

# Controles antes de continuar.
stopifnot(inherits(apt_climatica, "SpatRaster"))
stopifnot(terra::nlyr(apt_climatica) == 1)
# rango_aptitud <- terra::minmax(apt_climatica)
# print(rango_aptitud)
# Obtener los valores mínimo y máximo
rango_aptitud <- terra::minmax(apt_climatica)

# Construir una tabla organizada
tabla_rango_aptitud <- data.frame(
  Indicador = c(
    "Valor mínimo",
    "Valor máximo"
  ),
  `Índice de aptitud climática` = c(
    as.numeric(rango_aptitud[1, 1]),
    as.numeric(rango_aptitud[2, 1])
  ),
  Interpretación = c(
    "Ausencia de aptitud climática",
    "Aptitud climática óptima"
  ),
  check.names = FALSE
)


tabla_rango_aptitud <- normalizar_tabla_utf8(tabla_rango_aptitud)
# Presentar la tabla con formato profesional
knitr::kable(
  tabla_rango_aptitud,
  caption = "Rango del índice de aptitud climática para la caña de azúcar",
  col.names = c(
    "Indicador",
    "Índice de aptitud climática",
    "Interpretación"
  ),
  digits = 2,
  align = c("l", "c", "l")
) |>
  kableExtra::kable_styling(
    bootstrap_options = c(
      "striped",
      "hover",
      "condensed",
      "responsive"
    ),
    full_width = TRUE,
    position = "center"
  ) |>
  kableExtra::column_spec(
    column = 1,
    bold = TRUE,
    color = "#123F72"
  ) |>
  kableExtra::column_spec(
    column = 2,
    bold = TRUE,
    color = "white",
    background = "#1F5B96",
    width = "4cm"
  ) |>
  kableExtra::column_spec(
    column = 3,
    width = "7cm"
  )

Rango del índice de aptitud climática para la caña de azúcar

Indicador	Índice de aptitud climática	Interpretación
Valor mínimo	0	Ausencia de aptitud climática
Valor máximo	1	Aptitud climática óptima

clase_aptitud <- terra::classify(
  apt_climatica,
  rcl = matrix(c(
    -Inf, 0.25, 1,
    0.25, 0.50, 2,
    0.50, 0.75, 3,
    0.75, Inf, 4
  ), ncol = 3, byrow = TRUE),
  include.lowest = TRUE,
  right = FALSE
)

levels(clase_aptitud) <- data.frame(
  value = 1:4,
  clase = c("Muy baja", "Baja", "Media", "Alta")
)

terra::writeRaster(
  c(tmed_anual, prec_anual, apt_temp, apt_prec, apt_climatica),
  file.path(ruta_resultados, "aptitud_climatica_cana.tif"),
  overwrite = TRUE
)

3.1 Mapas globales

plot(apt_climatica,
     col = hcl.colors(100, "YlGn", rev = FALSE),
     range = c(0, 1),
     main = "Aptitud climática global para caña de azúcar",
     axes = TRUE)
lines(paises, col = "grey35", lwd = 0.25)

plot(clase_aptitud,
     col = c("#d73027", "#fee08b", "#91cf60", "#1a9850"),
     main = "Clases de aptitud climática")
lines(paises, col = "grey30", lwd = 0.25)

3.2 Países con mayor superficie de aptitud alta

En lugar de escoger países subjetivamente, se calcula la superficie de celdas con índice mayor o igual a 0,75. La superficie se estima con cellSize() en km² y se resume por país.

area_alta <- ifel(apt_climatica >= 0.75, cellSize(apt_climatica, unit = "km"), NA)
area_por_pais <- terra::extract(area_alta, paises, fun = sum, na.rm = TRUE)

campo_nombre <- intersect(c("NAME_0", "COUNTRY", "NAME_ENGLI", "ADMIN"), names(paises))[1]
if (is.na(campo_nombre)) stop("No se encontró el campo de nombre de país en la capa mundial.")

ranking_paises <- data.frame(
  pais = paises[[campo_nombre]],
  area_alta_km2 = area_por_pais[, 2]
) |>
  dplyr::filter(is.finite(area_alta_km2), area_alta_km2 > 0) |>
  dplyr::arrange(dplyr::desc(area_alta_km2)) |>
  dplyr::mutate(area_alta_km2 = round(area_alta_km2, 0))

knitr::kable(head(ranking_paises, 10),
             caption = "Diez países con mayor superficie climáticamente apta") |>
  kableExtra::kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                            full_width = TRUE)

Diez países con mayor superficie climáticamente apta

NAME_0	area_alta_km2
Brazil	5761923
Democratic Republic of the Congo	1800743
Indonesia	1602374
Colombia	860486
India	681546
Peru	668296
Venezuela	645191
Bolivia	462795
Central African Republic	392362
Cameroon	366904

paises_top3_nombres <- head(ranking_paises$pais, 3)
# El objeto se reconstruye después de aplicar, si fuera necesario, el respaldo
# por aptitud máxima en el bloque de mapas.
paises_top3 <- paises[paises[[campo_nombre]] %in% paises_top3_nombres, ]

3.3 Mapas de países con mayor superficie climática apta (Top 3)

# ==============================================================
# MAPAS DE PAÍSES CON MAYOR SUPERFICIE CLIMÁTICAMENTE APTA
# ==============================================================

# 1. Validar los objetos requeridos
if (!exists("apt_climatica")) {
  stop("No existe el objeto 'apt_climatica'.")
}

if (!exists("paises")) {
  stop("No existe el objeto espacial 'paises'.")
}

if (!inherits(apt_climatica, "SpatRaster")) {
  stop("'apt_climatica' no es un objeto SpatRaster.")
}

if (!inherits(paises, "SpatVector")) {
  stop("'paises' no es un objeto SpatVector.")
}

# 2. Extraer los atributos de la capa espacial
atributos_paises <- terra::values(
  paises,
  dataframe = TRUE
)

if (
  !is.data.frame(atributos_paises) ||
  nrow(atributos_paises) == 0L
) {
  stop(
    "La capa mundial no contiene una tabla de atributos válida."
  )
}

# 3. Identificar automáticamente la columna con el nombre del país
campos_posibles <- c(
  "NAME_0",
  "COUNTRY",
  "ADMIN",
  "NAME_ENGLI",
  "NAME",
  "SOVEREIGNT"
)

campos_encontrados <- intersect(
  campos_posibles,
  names(atributos_paises)
)

if (length(campos_encontrados) > 0L) {

  campo_nombre <- campos_encontrados[1]

} else {

  columnas_texto <- names(atributos_paises)[
    vapply(
      atributos_paises,
      function(x) {
        is.character(x) || is.factor(x)
      },
      logical(1)
    )
  ]

  if (length(columnas_texto) == 0L) {
    stop(
      paste0(
        "No se identificó una columna de nombres. ",
        "Columnas disponibles: ",
        paste(names(atributos_paises), collapse = ", ")
      )
    )
  }

  cantidad_unicos <- vapply(
    atributos_paises[columnas_texto],
    function(x) {
      length(
        unique(
          stats::na.omit(
            as.character(x)
          )
        )
      )
    },
    integer(1)
  )

  campo_nombre <- columnas_texto[
    which.max(cantidad_unicos)
  ]
}

# 4. Extraer los nombres de países
nombres_paises <- trimws(
  as.character(
    atributos_paises[[campo_nombre]]
  )
)

numero_paises <- length(nombres_paises)

if (numero_paises == 0L) {
  stop(
    paste0(
      "La columna '",
      campo_nombre,
      "' no contiene nombres de países."
    )
  )
}

if (numero_paises != terra::nrow(paises)) {
  stop(
    paste0(
      "La cantidad de nombres (",
      numero_paises,
      ") no coincide con la cantidad de geometrías (",
      terra::nrow(paises),
      ")."
    )
  )
}

# 5. Construir el catálogo, conservando el índice espacial
catalogo_paises <- data.frame(
  ID = seq_along(nombres_paises),
  pais = nombres_paises,
  stringsAsFactors = FALSE
)

# 6. Calcular el área de cada celda en km²
area_celda_km2 <- terra::cellSize(
  apt_climatica,
  unit = "km"
)

# 7. Conservar el área con aptitud climática alta
area_alta <- terra::ifel(
  apt_climatica >= 0.75,
  area_celda_km2,
  NA
)

names(area_alta) <- "area_alta_km2"

# 8. Calcular cada estadística por separado
# No se utilizan funciones personalizadas multivariadas.

resultado_area <- terra::extract(
  area_alta,
  paises,
  fun = "sum",
  na.rm = TRUE
)

resultado_media <- terra::extract(
  apt_climatica,
  paises,
  fun = "mean",
  na.rm = TRUE
)

resultado_maxima <- terra::extract(
  apt_climatica,
  paises,
  fun = "max",
  na.rm = TRUE
)

# 9. Verificar y normalizar las tablas extraídas
if (ncol(resultado_area) < 2L) {
  stop("No se pudo calcular el área altamente apta por país.")
}

if (ncol(resultado_media) < 2L) {
  stop("No se pudo calcular la aptitud media por país.")
}

if (ncol(resultado_maxima) < 2L) {
  stop("No se pudo calcular la aptitud máxima por país.")
}

resultado_area <- data.frame(
  ID = as.integer(resultado_area[[1]]),
  area_alta_km2 = as.numeric(resultado_area[[2]])
)

resultado_media <- data.frame(
  ID = as.integer(resultado_media[[1]]),
  aptitud_media = as.numeric(resultado_media[[2]])
)

resultado_maxima <- data.frame(
  ID = as.integer(resultado_maxima[[1]]),
  aptitud_maxima = as.numeric(resultado_maxima[[2]])
)

# 10. Integrar los resultados
estadisticas_paises <- catalogo_paises |>
  dplyr::left_join(
    resultado_area,
    by = "ID"
  ) |>
  dplyr::left_join(
    resultado_media,
    by = "ID"
  ) |>
  dplyr::left_join(
    resultado_maxima,
    by = "ID"
  )

# Los países sin área de aptitud alta reciben cero
estadisticas_paises$area_alta_km2[
  is.na(estadisticas_paises$area_alta_km2)
] <- 0

# 11. Validación del resultado integrado
stopifnot(
  nrow(estadisticas_paises) == numero_paises,
  is.numeric(estadisticas_paises$area_alta_km2),
  is.numeric(estadisticas_paises$aptitud_media),
  is.numeric(estadisticas_paises$aptitud_maxima)
)

# 12. Conservar países con nombres y datos climáticos válidos
estadisticas_validas <- estadisticas_paises |>
  dplyr::filter(
    !is.na(pais),
    nzchar(trimws(pais)),
    is.finite(aptitud_maxima)
  ) |>
  dplyr::arrange(
    dplyr::desc(area_alta_km2),
    dplyr::desc(aptitud_media),
    dplyr::desc(aptitud_maxima)
  )

if (nrow(estadisticas_validas) == 0L) {
  stop(
    "No existen países con estadísticas climáticas válidas."
  )
}

# 13. Seleccionar los tres primeros países
top_paises <- estadisticas_validas |>
  dplyr::slice_head(
    n = min(3L, nrow(estadisticas_validas))
  )

paises_top3_nombres <- top_paises$pais
n_mapas <- nrow(top_paises)

# cat(
#  "Campo utilizado para el nombre del país:",
#  campo_nombre,
#  "\n"
# )

# cat(
#  "Países seleccionados:",
#  paste(paises_top3_nombres, collapse = ", "),
#  "\n"
# )

# Traducir los nombres para la presentación del informe
traduccion_paises <- c(
  "Brazil" = "Brasil",
  "Democratic Republic of the Congo" = "República Democrática del Congo",
  "Indonesia" = "Indonesia"
)

paises_presentacion <- ifelse(
  paises_top3_nombres %in% names(traduccion_paises),
  traduccion_paises[paises_top3_nombres],
  paises_top3_nombres
)

# Construir un recuadro HTML profesional
html_paises <- paste0(
  '<div class="recuadro-seleccion">',
  '<div class="recuadro-titulo">',
  'Países seleccionados para el análisis detallado',
  '</div>',
  '<div class="recuadro-contenido">',
  paste(
    paste0(
      '<span class="etiqueta-pais">',
      paises_presentacion,
      '</span>'
    ),
    collapse = " "
  ),
  '</div>',
  '<div class="recuadro-nota">',
  'Selección basada en la superficie estimada con aptitud climática alta.',
  '</div>',
  '</div>'
)

htmltools::HTML(html_paises)

Países seleccionados para el análisis detallado

Brasil República Democrática del Congo Indonesia

Selección basada en la superficie estimada con aptitud climática alta.

# 14. Tabla de resultados
tabla_top_paises <- top_paises |>
  dplyr::transmute(
    País = pais,
    `Área con aptitud alta (km²)` = round(
      area_alta_km2,
      0
    ),
    `Aptitud media` = round(
      aptitud_media,
      3
    ),
    `Aptitud máxima` = round(
      aptitud_maxima,
      3
    )
  )

tabla_top_paises<- normalizar_tabla_utf8(tabla_top_paises)
tabla_resultados <- knitr::kable(
  tabla_top_paises,
  caption = paste(
    "Países con mayor superficie climáticamente",
    "apta para caña de azúcar"
  ),
  align = c("l", "r", "c", "c"),
  format.args = list(
    big.mark = ".",
    decimal.mark = ","
  )
) |>
  kableExtra::kable_styling(
    bootstrap_options = c(
      "striped",
      "hover",
      "condensed",
      "responsive"
    ),
    full_width = TRUE,
    position = "center"
  )

# print(tabla_resultados)
tabla_resultados

Países con mayor superficie climáticamente apta para caña de azúcar

País	Área con aptitud alta (km²)	Aptitud media	Aptitud máxima
Brazil	5.761.923	0,756	1
Democratic Republic of the Congo	1.800.743	0,860	1
Indonesia	1.602.374	0,904	1

# 15. Guardar la configuración gráfica
configuracion_anterior <- graphics::par(
  no.readonly = TRUE
)

# 16. Definir la disposición de los mapas
graphics::par(
  mfrow = c(n_mapas, 1),
  mar = c(3.5, 3.5, 4.5, 5.5),
  oma = c(0, 0, 2, 0)
)

# 17. Crear los mapas
for (i in seq_len(n_mapas)) {

  id_pais <- as.integer(
    top_paises$ID[i]
  )

  nombre_pais <- as.character(
    top_paises$pais[i]
  )

  if (
    !is.finite(id_pais) ||
    id_pais < 1L ||
    id_pais > terra::nrow(paises)
  ) {
    graphics::plot.new()

    graphics::text(
      0.5,
      0.5,
      paste(
        "Identificador espacial inválido para",
        nombre_pais
      )
    )

  } else {

    poligono_pais <- paises[id_pais, ]

    raster_pais <- terra::crop(
      apt_climatica,
      poligono_pais,
      snap = "out"
    )

    raster_pais <- terra::mask(
      raster_pais,
      poligono_pais
    )

    cantidad_celdas <- terra::global(
      !is.na(raster_pais),
      fun = "sum",
      na.rm = TRUE
    )[1, 1]

    if (
      is.na(cantidad_celdas) ||
      cantidad_celdas == 0
    ) {
      graphics::plot.new()

      graphics::text(
        0.5,
        0.5,
        paste(
          "No existen datos climáticos válidos para",
          nombre_pais
        )
      )

    } else {

      terra::plot(
        raster_pais,
        col = hcl.colors(
          100,
          "YlGn"
        ),
        range = c(0, 1),
        axes = TRUE,
        main = paste(
          "Aptitud climática para caña de azúcar:",
          nombre_pais
        ),
        plg = list(
          title = "Índice de aptitud",
          cex = 0.8
        )
      )

      terra::lines(
        poligono_pais,
        col = "#263238",
        lwd = 1.1
      )
    }
  }
}

# 18. Título general
graphics::mtext(
  text = "Países con mayor superficie climáticamente apta",
  side = 3,
  outer = TRUE,
  line = 0.3,
  cex = 1.15,
  font = 2,
  col = "#123f72"
)

# 19. Restaurar la configuración gráfica
graphics::par(
  configuracion_anterior
)

Interpretación. Los tres países anteriores no se fijan de antemano: son el resultado del área clasificada como alta bajo la resolución, periodo y umbrales seleccionados. La presencia de aptitud climática no implica que toda esa superficie esté disponible para agricultura; pueden existir bosques, áreas protegidas, ciudades, pendientes fuertes o suelos no aptos.

4 Sitios de referencia en el Valle del Cauca

Se usan tres puntos representativos del corredor geográfico del valle: Palmira, Candelaria y Zarzal. Son puntos urbanos de referencia, no parcelas agrícolas. Las ubicaciones se representan directamente en un mapa interactivo elaborado dentro del informe. El mapa utiliza teselas de OpenStreetMap, no requiere claves API y muestra los tres sitios con etiquetas permanentes.

sitios <- data.frame(
  sitio = c("Palmira", "Candelaria", "Zarzal"),
  longitud = c(-76.3036, -76.3482, -76.0715),
  latitud  = c(  3.5394,   3.4067,   4.3946)
)

sitios_v <- vect(sitios, geom = c("longitud", "latitud"), crs = "EPSG:4326")
sitios <- normalizar_tabla_utf8(sitios)
knitr::kable(sitios, caption = "Sitios climáticos de referencia en el Valle del Cauca") |>
  kableExtra::kable_styling(bootstrap_options = c("striped", "hover"), full_width = TRUE)

Sitios climáticos de referencia en el Valle del Cauca

sitio	longitud	latitud
Palmira	-76.3036	3.5394
Candelaria	-76.3482	3.4067
Zarzal	-76.0715	4.3946

# Mapa interactivo sin dependencia de CARTO ni claves API.
# OpenStreetMap se agrega mediante su URL oficial y atribución visible.
paleta_sitios <- leaflet::colorFactor(
  palette = c("#123f72", "#c7a34b", "#2f855a"),
  domain = sitios$sitio
)

mapa_valle <- leaflet::leaflet(
  data = sitios,
  width = "100%",
  height = 560,
  options = leaflet::leafletOptions(minZoom = 7, maxZoom = 18)
) |>
  leaflet::addTiles(
    urlTemplate = "https://tile.openstreetmap.org/{z}/{x}/{y}.png",
    attribution = "&copy; OpenStreetMap contributors",
    options = leaflet::tileOptions(maxZoom = 19)
  ) |>
  leaflet::addCircleMarkers(
    lng = ~longitud,
    lat = ~latitud,
    radius = 9,
    stroke = TRUE,
    weight = 3,
    color = "#ffffff",
    fillColor = ~paleta_sitios(sitio),
    fillOpacity = 1,
    label = ~sitio,
    labelOptions = leaflet::labelOptions(
      noHide = TRUE,
      direction = "top",
      offset = c(0, -10),
      textOnly = FALSE,
      style = list(
        "font-weight" = "bold",
        "color" = "#123f72",
        "background-color" = "white",
        "border" = "1px solid #123f72",
        "padding" = "3px 6px"
      )
    ),
    popup = ~paste0(
      "<strong>", sitio, "</strong><br>",
      "Latitud: ", sprintf("%.4f", latitud), "<br>",
      "Longitud: ", sprintf("%.4f", longitud)
    )
  ) |>
  leaflet::fitBounds(
    lng1 = min(sitios$longitud) - 0.12,
    lat1 = min(sitios$latitud) - 0.12,
    lng2 = max(sitios$longitud) + 0.12,
    lat2 = max(sitios$latitud) + 0.12
  ) |>
  leaflet::addScaleBar(position = "bottomleft")

mapa_valle

4.1 Extracción climática y series mensuales

WorldClim es una climatología mensual 1970–2000. Por tanto, las gráficas muestran el ciclo promedio enero–diciembre, no una serie interanual de observaciones.

clima_mensual <- c(tmed_m, prec_m)
ext_sitios <- terra::extract(clima_mensual, sitios_v)

serie_temp <- ext_sitios[, 2:13] |>
  dplyr::mutate(sitio = sitios$sitio) |>
  tidyr::pivot_longer(-sitio, names_to = "variable", values_to = "temperatura_C") |>
  dplyr::mutate(mes = rep(1:12, times = nrow(sitios)))

serie_prec <- ext_sitios[, 14:25] |>
  dplyr::mutate(sitio = sitios$sitio) |>
  tidyr::pivot_longer(-sitio, names_to = "variable", values_to = "precipitacion_mm") |>
  dplyr::mutate(mes = rep(1:12, times = nrow(sitios)))

serie_sitios <- dplyr::left_join(
  dplyr::select(serie_temp, sitio, mes, temperatura_C),
  dplyr::select(serie_prec, sitio, mes, precipitacion_mm),
  by = c("sitio", "mes")
) |>
  dplyr::mutate(mes_nombre = factor(month.abb[mes], levels = month.abb))

serie_sitios <- normalizar_tabla_utf8(serie_sitios)
knitr::kable(serie_sitios, digits = 1,
             caption = "Climatología mensual extraída en los sitios") |>
  kableExtra::kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
                            full_width = TRUE)

Climatología mensual extraída en los sitios

sitio	mes	temperatura_C	precipitacion_mm	mes_nombre
Palmira	1	22.5	92	Jan
Palmira	2	22.6	102	Feb
Palmira	3	22.7	132	Mar
Palmira	4	22.4	172	Apr
Palmira	5	22.4	146	May
Palmira	6	22.3	112	Jun
Palmira	7	22.6	61	Jul
Palmira	8	22.7	73	Aug
Palmira	9	22.6	105	Sep
Palmira	10	22.0	195	Oct
Palmira	11	21.9	173	Nov
Palmira	12	22.1	121	Dec
Candelaria	1	23.8	95	Jan
Candelaria	2	23.9	110	Feb
Candelaria	3	24.0	142	Mar
Candelaria	4	23.8	185	Apr
Candelaria	5	23.7	168	May
Candelaria	6	23.7	106	Jun
Candelaria	7	24.0	55	Jul
Candelaria	8	24.1	78	Aug
Candelaria	9	23.9	112	Sep
Candelaria	10	23.4	210	Oct
Candelaria	11	23.3	176	Nov
Candelaria	12	23.5	128	Dec
Zarzal	1	23.5	81	Jan
Zarzal	2	23.7	102	Feb
Zarzal	3	23.8	141	Mar
Zarzal	4	23.5	193	Apr
Zarzal	5	23.4	182	May
Zarzal	6	23.3	127	Jun
Zarzal	7	23.8	74	Jul
Zarzal	8	23.8	105	Aug
Zarzal	9	23.5	135	Sep
Zarzal	10	22.8	214	Oct
Zarzal	11	22.7	168	Nov
Zarzal	12	23.0	110	Dec

g_temp <- ggplot(serie_sitios, aes(mes_nombre, temperatura_C, color = sitio, group = sitio)) +
  geom_line(linewidth = 0.9) + geom_point(size = 1.8) +
  scale_color_viridis_d(option = "C", end = 0.85) +
  labs(title = "Ciclo climatológico mensual de temperatura",
       x = "Mes", y = "Temperatura media (°C)", color = "Sitio") +
  theme_minimal(base_size = 11) + theme(legend.position = "bottom")

g_prec <- ggplot(serie_sitios, aes(mes_nombre, precipitacion_mm, color = sitio, group = sitio)) +
  geom_line(linewidth = 0.9) + geom_point(size = 1.8) +
  scale_color_viridis_d(option = "D", end = 0.85) +
  labs(title = "Ciclo climatológico mensual de precipitación",
       x = "Mes", y = "Precipitación (mm/mes)", color = "Sitio") +
  theme_minimal(base_size = 11) + theme(legend.position = "bottom")

g_temp / g_prec

5 Similaridad climática global

La distancia euclidiana directa entre variables con unidades diferentes estaría dominada por la precipitación. Para evitarlo, las 24 variables mensuales se estandarizan con media 0 y desviación estándar 1. Para cada sitio se calcula:

\[ d_i = \sqrt{\sum_{j=1}^{24}(z_{ij}-z_{sj})^2} \]

Luego se transforma la distancia en una similaridad acotada entre 0 y 1:

\[ S_i = \exp(-d_i) \]

Valores cercanos a 1 representan mayor analogía climática con el sitio del Valle del Cauca.

clima_24 <- c(tmed_m, prec_m)

# Estadísticos globales por capa; se aplican los mismos parámetros a raster y sitios.
medias <- as.numeric(global(clima_24, "mean", na.rm = TRUE)[, 1])
desv <- as.numeric(global(clima_24, "sd", na.rm = TRUE)[, 1])
if (any(!is.finite(desv) | desv == 0)) stop("Hay capas sin variabilidad válida.")

clima_z <- (clima_24 - medias) / desv
sitios_24 <- as.matrix(ext_sitios[, 2:25])
sitios_z <- sweep(sweep(sitios_24, 2, medias, "-"), 2, desv, "/")

mapas_sim <- list()
for (i in seq_len(nrow(sitios))) {
  dif2 <- (clima_z - sitios_z[i, ])^2
  distancia <- sqrt(sum(dif2, na.rm = FALSE))
  similitud <- exp(-distancia)
  names(similitud) <- paste0("sim_", sitios$sitio[i])
  mapas_sim[[i]] <- similitud
}
sim_global <- rast(mapas_sim)
writeRaster(sim_global, file.path(ruta_resultados, "similaridad_sitios_valle.tif"), overwrite = TRUE)

par(mfrow = c(3, 1), mar = c(2, 3, 3, 5))
for (i in 1:3) {
  plot(sim_global[[i]], col = hcl.colors(100, "Inferno", rev = TRUE),
       range = c(0, 1), main = paste("Similaridad climática con", sitios$sitio[i]),
       axes = TRUE)
  lines(paises, col = "grey70", lwd = 0.2)
}

par(mfrow = c(1, 1))

5.1 Consenso entre sitios

El promedio resume la analogía general con los tres sitios; el mínimo es más exigente y conserva únicamente zonas semejantes a todos ellos.

sim_promedio <- mean(sim_global, na.rm = TRUE)
sim_conservadora <- terra::ifel(
  sim_global[[1]] <= sim_global[[2]] & sim_global[[1]] <= sim_global[[3]], sim_global[[1]],
  terra::ifel(sim_global[[2]] <= sim_global[[3]], sim_global[[2]], sim_global[[3]])
)
names(sim_promedio) <- "similitud_promedio"
names(sim_conservadora) <- "similitud_conservadora"

par(mfrow = c(2, 1), mar = c(2, 3, 3, 5))
plot(sim_promedio, col = hcl.colors(100, "Inferno", rev = TRUE), range = c(0, 1),
     main = "Similaridad promedio con los sitios del Valle del Cauca")
lines(paises, col = "grey70", lwd = 0.2)
plot(sim_conservadora, col = hcl.colors(100, "Inferno", rev = TRUE), range = c(0, 1),
     main = "Similaridad conservadora: mínimo entre los tres sitios")
lines(paises, col = "grey70", lwd = 0.2)

par(mfrow = c(1, 1))

# ==============================================================
# SIMILARIDAD CLIMÁTICA GLOBAL
# CHUNK AUTÓNOMO: NO DEPENDE DE OBJETOS PREVIOS EN MEMORIA
# ==============================================================

construir_similaridad_climatica <- function() {

  # ------------------------------------------------------------
  # 1. Verificar paquetes
  # ------------------------------------------------------------

  paquetes_necesarios <- c(
    "terra",
    "geodata",
    "knitr",
    "kableExtra"
  )

  paquetes_faltantes <- paquetes_necesarios[
    !vapply(
      paquetes_necesarios,
      requireNamespace,
      logical(1),
      quietly = TRUE
    )
  ]

  if (length(paquetes_faltantes) > 0L) {
    stop(
      paste0(
        "Faltan los paquetes: ",
        paste(paquetes_faltantes, collapse = ", "),
        ". Instálelos antes de renderizar."
      ),
      call. = FALSE
    )
  }

  terra::terraOptions(progress = 0)

  # ------------------------------------------------------------
  # 2. Crear las carpetas necesarias
  # ------------------------------------------------------------

  ruta_datos_local <- file.path(
    getwd(),
    "datos_clima"
  )

  ruta_resultados_local <- file.path(
    getwd(),
    "resultados_informe2"
  )

  dir.create(
    ruta_datos_local,
    recursive = TRUE,
    showWarnings = FALSE
  )

  dir.create(
    ruta_resultados_local,
    recursive = TRUE,
    showWarnings = FALSE
  )

  # ------------------------------------------------------------
  # 3. Cargar o descargar WorldClim
  # ------------------------------------------------------------

  tmed_local <- geodata::worldclim_global(
    var = "tavg",
    res = 10,
    path = ruta_datos_local
  )

  prec_local <- geodata::worldclim_global(
    var = "prec",
    res = 10,
    path = ruta_datos_local
  )

  if (!inherits(tmed_local, "SpatRaster")) {
    stop(
      "No fue posible cargar la temperatura mensual.",
      call. = FALSE
    )
  }

  if (!inherits(prec_local, "SpatRaster")) {
    stop(
      "No fue posible cargar la precipitación mensual.",
      call. = FALSE
    )
  }

  if (terra::nlyr(tmed_local) != 12L) {
    stop(
      paste0(
        "La temperatura debe contener 12 capas, pero contiene ",
        terra::nlyr(tmed_local),
        "."
      ),
      call. = FALSE
    )
  }

  if (terra::nlyr(prec_local) != 12L) {
    stop(
      paste0(
        "La precipitación debe contener 12 capas, pero contiene ",
        terra::nlyr(prec_local),
        "."
      ),
      call. = FALSE
    )
  }

  if (
    !isTRUE(
      terra::compareGeom(
        tmed_local,
        prec_local,
        stopOnError = FALSE
      )
    )
  ) {
    stop(
      paste(
        "Las capas de temperatura y precipitación",
        "no son espacialmente compatibles."
      ),
      call. = FALSE
    )
  }

  names(tmed_local) <- paste0(
    "tmed_",
    sprintf("%02d", 1:12)
  )

  names(prec_local) <- paste0(
    "prec_",
    sprintf("%02d", 1:12)
  )

  # WorldClim puede almacenar temperatura con factor 10.
  temperatura_prueba <- terra::global(
    tmed_local[[1]],
    fun = "mean",
    na.rm = TRUE
  )[1, 1]

  if (
    is.finite(temperatura_prueba) &&
    temperatura_prueba > 60
  ) {
    tmed_local <- tmed_local / 10
  }

  # ------------------------------------------------------------
  # 4. Construir las 24 capas climáticas
  # ------------------------------------------------------------

  clima_local <- c(
    tmed_local,
    prec_local
  )

  if (terra::nlyr(clima_local) != 24L) {
    stop(
      "No fue posible construir las 24 capas climáticas.",
      call. = FALSE
    )
  }

  # ------------------------------------------------------------
  # 5. Crear autónomamente los sitios del Valle del Cauca
  # ------------------------------------------------------------

  sitios_local <- data.frame(
    sitio = c(
      "Palmira",
      "Candelaria",
      "Zarzal"
    ),
    longitud = c(
      -76.3036,
      -76.3482,
      -76.0715
    ),
    latitud = c(
      3.5394,
      3.4067,
      4.3946
    ),
    stringsAsFactors = FALSE
  )

  sitios_vector <- terra::vect(
    sitios_local,
    geom = c(
      "longitud",
      "latitud"
    ),
    crs = "EPSG:4326"
  )

  # ------------------------------------------------------------
  # 6. Extraer el clima de los sitios
  # ------------------------------------------------------------

  extraccion_local <- terra::extract(
    clima_local,
    sitios_vector
  )

  if (
    !is.data.frame(extraccion_local) ||
    nrow(extraccion_local) != nrow(sitios_local) ||
    ncol(extraccion_local) != 25L
  ) {
    stop(
      paste0(
        "La extracción climática no tiene la estructura esperada. ",
        "Dimensiones obtenidas: ",
        nrow(extraccion_local),
        " × ",
        ncol(extraccion_local),
        "."
      ),
      call. = FALSE
    )
  }

  sitios_clima <- as.matrix(
    extraccion_local[
      ,
      2:25,
      drop = FALSE
    ]
  )

  storage.mode(sitios_clima) <- "numeric"

  if (any(!is.finite(sitios_clima))) {
    stop(
      "La extracción climática de los sitios contiene valores no válidos.",
      call. = FALSE
    )
  }

  # ------------------------------------------------------------
  # 7. Calcular estadísticas globales por capa
  # ------------------------------------------------------------

  medias_tabla <- terra::global(
    clima_local,
    fun = "mean",
    na.rm = TRUE
  )

  desviaciones_tabla <- terra::global(
    clima_local,
    fun = "sd",
    na.rm = TRUE
  )

  medias_locales <- as.numeric(
    medias_tabla[, 1]
  )

  desviaciones_locales <- as.numeric(
    desviaciones_tabla[, 1]
  )

  if (
    length(medias_locales) != 24L ||
    length(desviaciones_locales) != 24L
  ) {
    stop(
      "No se obtuvieron 24 medias y 24 desviaciones estándar.",
      call. = FALSE
    )
  }

  if (any(!is.finite(medias_locales))) {
    stop(
      "Existen medias climáticas no válidas.",
      call. = FALSE
    )
  }

  if (
    any(
      !is.finite(desviaciones_locales) |
      desviaciones_locales <= 0
    )
  ) {
    stop(
      "Existen capas climáticas sin variabilidad válida.",
      call. = FALSE
    )
  }

  # `terra::global()` devuelve una fila por capa.
  # Las funciones incorporadas mean y sd son adecuadas
  # para rásteres globales grandes. 

  # ------------------------------------------------------------
  # 8. Estandarizar las capas climáticas
  # ------------------------------------------------------------

  capas_estandarizadas <- vector(
    mode = "list",
    length = 24L
  )

  for (j in seq_len(24L)) {

    capas_estandarizadas[[j]] <- (
      clima_local[[j]] - medias_locales[j]
    ) / desviaciones_locales[j]
  }

  clima_estandarizado <- terra::rast(
    capas_estandarizadas
  )

  names(clima_estandarizado) <- paste0(
    "z_",
    names(clima_local)
  )

  if (terra::nlyr(clima_estandarizado) != 24L) {
    stop(
      "No fue posible estandarizar las 24 capas.",
      call. = FALSE
    )
  }

  # ------------------------------------------------------------
  # 9. Estandarizar el clima de los sitios
  # ------------------------------------------------------------

  sitios_estandarizados <- sweep(
    sitios_clima,
    MARGIN = 2,
    STATS = medias_locales,
    FUN = "-"
  )

  sitios_estandarizados <- sweep(
    sitios_estandarizados,
    MARGIN = 2,
    STATS = desviaciones_locales,
    FUN = "/"
  )

  if (any(!is.finite(sitios_estandarizados))) {
    stop(
      "La estandarización climática de los sitios produjo valores no válidos.",
      call. = FALSE
    )
  }

  # ------------------------------------------------------------
  # 10. Calcular la distancia euclidiana
  # ------------------------------------------------------------

  mapas_locales <- vector(
    mode = "list",
    length = nrow(sitios_local)
  )

  for (i in seq_len(nrow(sitios_local))) {

    # Inicializar con la primera variable
    distancia_cuadrada <- (
      clima_estandarizado[[1]] -
      sitios_estandarizados[i, 1]
    )^2

    # Acumular las otras 23 diferencias cuadráticas
    for (j in 2:24) {

      distancia_cuadrada <- distancia_cuadrada + (
        clima_estandarizado[[j]] -
        sitios_estandarizados[i, j]
      )^2
    }

    distancia_euclidiana <- sqrt(
      distancia_cuadrada
    )

    mapa_similaridad <- exp(
      -distancia_euclidiana
    )

    names(mapa_similaridad) <- paste0(
      "sim_",
      make.names(sitios_local$sitio[i])
    )

    mapas_locales[[i]] <- mapa_similaridad
  }

  # ------------------------------------------------------------
  # 11. Integrar y validar los tres mapas
  # ------------------------------------------------------------

  if (
    length(mapas_locales) != 3L ||
    !all(
      vapply(
        mapas_locales,
        inherits,
        logical(1),
        what = "SpatRaster"
      )
    )
  ) {
    stop(
      "No fue posible construir los tres mapas de similaridad.",
      call. = FALSE
    )
  }

  similitud_global <- c(
    mapas_locales[[1]],
    mapas_locales[[2]],
    mapas_locales[[3]]
  )

  names(similitud_global) <- paste0(
    "Similaridad_",
    make.names(sitios_local$sitio)
  )

  if (terra::nlyr(similitud_global) != 3L) {
    stop(
      "El resultado final no contiene tres mapas.",
      call. = FALSE
    )
  }

  rango_local <- terra::minmax(
    similitud_global
  )

  if (
    any(rango_local < 0, na.rm = TRUE) ||
    any(rango_local > 1, na.rm = TRUE)
  ) {
    stop(
      "La similaridad contiene valores fuera del intervalo 0–1.",
      call. = FALSE
    )
  }

   # ------------------------------------------------------------
  # 12. Crear y mostrar la tabla de resumen
  # ------------------------------------------------------------

  resumen_local <- data.frame(
    Sitio = sitios_local$sitio,
    `Similaridad mínima` = round(
      as.numeric(rango_local[1, ]),
      4
    ),
    `Similaridad máxima` = round(
      as.numeric(rango_local[2, ]),
      4
    ),
    check.names = FALSE,
    stringsAsFactors = FALSE
  )
  
assign(
  "resumen_similaridad",
  resumen_local,
  envir = parent.frame()
)

#  tabla_local <- knitr::kable(
#    resumen_local,
#    format = "html",
#    caption = paste(
#      "Rango de similaridad climática global",
#      "por sitio del Valle del Cauca"
#    ),
#    col.names = c(
#      "Sitio del Valle del Cauca",
#      "Similaridad mínima",
#      "Similaridad máxima"
#    ),
#    align = c("l", "c", "c"),
#    digits = 4,
#    escape = TRUE,
#    table.attr = 'class="table table-striped table-hover table-condensed"'
#  ) |>
#    kableExtra::kable_styling(
#      bootstrap_options = c(
#        "striped",
#        "hover",
#        "condensed",
#        "responsive"
#      ),
#      full_width = TRUE,
#      position = "center",
#      font_size = 14
#    ) |>
#    kableExtra::row_spec(
#      row = 0,
#      bold = TRUE,
#      color = "white",
#      background = "#123F72"
#    ) |>
#    kableExtra::column_spec(
#      column = 1,
#      bold = TRUE,
#      color = "#123F72",
#      width = "35%"
#    ) |>
#    kableExtra::column_spec(
#      column = 2,
#      width = "25%"
#    ) |>
#    kableExtra::column_spec(
#      column = 3,
#      width = "25%"
#    )

  # Emitir el HTML como contenido del informe, no como texto de consola.
#  cat(
#    as.character(tabla_local),
#    sep = "\n"
#  )

  # Separación entre la tabla y los mapas.
#  cat(
#    "\n<div style='margin-bottom: 28px;'></div>\n"
# )
  
  # ------------------------------------------------------------
  # 13. Cargar límites mundiales
  # ------------------------------------------------------------

  paises_locales <- geodata::world(
    path = ruta_datos_local
  )

  # ------------------------------------------------------------
  # 14. Graficar los mapas
  # ------------------------------------------------------------

  configuracion_anterior <- graphics::par(
    no.readonly = TRUE
  )

  graphics::par(
    mfrow = c(3, 1),
    mar = c(2.5, 3, 4, 5.5),
    oma = c(0, 0, 2, 0)
  )

  for (i in seq_len(3L)) {

    terra::plot(
      similitud_global[[i]],
      col = hcl.colors(
        100,
        palette = "Inferno",
        rev = TRUE
      ),
      range = c(0, 1),
      axes = TRUE,
      main = paste(
        "Similaridad climática global con",
        sitios_local$sitio[i]
      ),
      plg = list(
        title = "Similaridad",
        cex = 0.8
      )
    )

    terra::lines(
      paises_locales,
      col = "grey55",
      lwd = 0.25
    )
  }

  graphics::mtext(
    text = paste(
      "Analogía climática con sitios",
      "del Valle del Cauca"
    ),
    side = 3,
    outer = TRUE,
    line = 0.3,
    font = 2,
    cex = 1.15,
    col = "#123F72"
  )

  graphics::par(
    configuracion_anterior
  )

  # ------------------------------------------------------------
  # 15. Guardar resultados
  # ------------------------------------------------------------

  archivo_similaridad <- file.path(
    ruta_resultados_local,
    "similaridad_climatica_sitios_valle.tif"
  )

  terra::writeRaster(
    similitud_global,
    filename = archivo_similaridad,
    overwrite = TRUE
  )

  # ------------------------------------------------------------
  # 16. Crear objetos requeridos por chunks posteriores
  # ------------------------------------------------------------

  entorno_salida <- parent.frame()

  assign(
    "tmed_m",
    tmed_local,
    envir = entorno_salida
  )

  assign(
    "prec_m",
    prec_local,
    envir = entorno_salida
  )

  assign(
    "clima_24",
    clima_local,
    envir = entorno_salida
  )

  assign(
    "clima_z",
    clima_estandarizado,
    envir = entorno_salida
  )

  assign(
    "sitios",
    sitios_local,
    envir = entorno_salida
  )

  assign(
    "ext_sitios",
    extraccion_local,
    envir = entorno_salida
  )

  assign(
    "sitios_z",
    sitios_estandarizados,
    envir = entorno_salida
  )

  assign(
    "paises",
    paises_locales,
    envir = entorno_salida
  )

  assign(
    "sim_global",
    similitud_global,
    envir = entorno_salida
  )

  assign(
    "mapas_sim",
    mapas_locales,
    envir = entorno_salida
  )

  assign(
    "ruta_resultados",
    ruta_resultados_local,
    envir = entorno_salida
  )

  invisible(
    similitud_global
  )
}

# ÚNICA LLAMADA DEL CHUNK
# No ejecutar las líneas internas de la función por separado.
construir_similaridad_climatica()

Rango de similaridad climática global por sitio del Valle del Cauca

Sitio del Valle del Cauca	Similaridad mínima	Similaridad máxima
Palmira	0	1
Candelaria	0	1
Zarzal	0	1

6 Comparación de aproximaciones

La aproximación por umbrales responde: ¿dónde se cumplen simultáneamente los rangos climáticos considerados óptimos? La analogía responde: ¿dónde se parece el ciclo mensual al de los sitios de referencia? No son equivalentes: dos lugares pueden compartir promedios anuales pero diferir en la estacionalidad.

# Escala relativa de similaridad alta: percentil 90 sobre celdas terrestres.
umbral_sim <- as.numeric(global(sim_promedio, function(x) quantile(x, 0.90, na.rm = TRUE))[, 1])

# Ráster numérico de clases. Se conserva sin niveles para los cálculos.
acuerdo_num <- terra::ifel(
  apt_climatica >= 0.75 & sim_promedio >= umbral_sim, 3,
  terra::ifel(
    apt_climatica >= 0.75, 1,
    terra::ifel(sim_promedio >= umbral_sim, 2, 0)
  )
)
names(acuerdo_num) <- "codigo_acuerdo"

# Se crea una copia categórica exclusivamente para la representación gráfica.
acuerdo_plot <- acuerdo_num
levels(acuerdo_plot) <- data.frame(
  value = 0:3,
  clase = c("Ninguna", "Solo umbrales", "Solo analogía", "Coincidencia")
)

plot(
  acuerdo_plot,
  col = c("grey90", "#2c7bb6", "#fdae61", "#1a9641"),
  main = "Comparación espacial de las dos aproximaciones"
)
lines(paises, col = "grey30", lwd = 0.2)

# Resumen robusto de superficies. Se evita zonal() sobre un ráster categórico,
# porque su primera columna puede devolverse como texto según la versión de terra.
area_celda <- terra::cellSize(acuerdo_num, unit = "km")
clases_codigo <- 0:3
clases_nombre <- c("Ninguna", "Solo umbrales", "Solo analogía", "Coincidencia")

area_por_clase <- vapply(
  clases_codigo,
  function(codigo) {
    area_clase <- terra::ifel(acuerdo_num == codigo, area_celda, NA)
    valor <- terra::global(area_clase, "sum", na.rm = TRUE)[1, 1]
    if (is.na(valor)) 0 else as.numeric(valor)
  },
  numeric(1)
)

resumen_acuerdo <- data.frame(
  codigo = clases_codigo,
  clase = clases_nombre,
  area_km2 = area_por_clase
)

stopifnot(
  nrow(resumen_acuerdo) == 4L,
  is.numeric(resumen_acuerdo$area_km2),
  all(is.finite(resumen_acuerdo$area_km2)),
  all(resumen_acuerdo$area_km2 >= 0)
)

resumen_acuerdo <- normalizar_tabla_utf8(resumen_acuerdo)
knitr::kable(
  resumen_acuerdo[, c("clase", "area_km2")],
  digits = 0,
  col.names = c("Clase", "Área aproximada (km²)"),
  caption = "Superficie por tipo de acuerdo"
)

Superficie por tipo de acuerdo

Clase	Área aproximada (km²)
Ninguna	115851391
Solo umbrales	14614790
Solo analogía	16788328
Coincidencia	7094662

6.1 Interpretación

• Coincidencia: zonas prioritarias para una evaluación agroecológica posterior, porque cumplen los rangos anuales y reproducen relativamente bien el ciclo mensual de los sitios vallecaucanos.
• Solo umbrales: el promedio anual es favorable, pero la distribución intraanual de lluvia o temperatura difiere del patrón de referencia.
• Solo analogía: el patrón mensual se parece, aunque al menos uno de los valores anuales no alcanza el criterio estricto de aptitud alta.
• Ninguna: baja prioridad climática bajo los supuestos adoptados.

7 Conclusiones

1. El análisis por umbrales permite localizar de manera transparente las celdas donde la temperatura media y la precipitación anual se aproximan a los rangos óptimos definidos para el caso.
2. El ranking de países se basa en superficie apta calculada, no en una selección visual. Los tres países y sus mapas se generan automáticamente al ejecutar el informe.
3. La analogía climática añade información sobre estacionalidad. Por ello puede excluir zonas con promedios anuales adecuados, pero con calendarios de lluvia o temperatura muy diferentes a los del Valle del Cauca.
4. Las áreas donde coinciden ambas aproximaciones constituyen candidatos razonables para estudios posteriores, pero no prueban rendimiento potencial en kg/ha.
5. Antes de una decisión productiva deben incorporarse suelo, pendiente, disponibilidad de agua, radiación, riesgo de inundación o sequía, infraestructura, restricciones ambientales, variedad y manejo. La etapa siguiente debe usar un modelo de cultivo calibrado y datos locales de rendimiento.

8 Limitaciones y reproducibilidad

• WorldClim 2.1 representa promedios de 1970–2000 y no las condiciones observadas de un año específico.
• La resolución de 10 minutos suaviza variaciones locales, especialmente en relieves complejos.
• Los tres puntos son referencias espaciales verificables mediante Google Maps, pero no sustituyen coordenadas de lotes o estaciones agrometeorológicas.
• Los rangos óptimos fueron fijados por el enunciado. Los márgenes externos del índice trapezoidal son supuestos metodológicos explícitos y pueden ajustarse.
• La similaridad depende de las variables, estandarización, métrica y transformación elegidas.
• La selección del percentil 90 se usa solo para comparar mapas; no representa un umbral agronómico universal.

write.csv(ranking_paises,
          file.path(ruta_resultados, "ranking_paises_aptitud.csv"), row.names = FALSE)
write.csv(serie_sitios,
          file.path(ruta_resultados, "clima_mensual_sitios_valle.csv"), row.names = FALSE)
cat("Resultados guardados en:", normalizePath(ruta_resultados), "\n")

Resultados guardados en: C:\Users\dcruz\OneDrive - Consorcio Programa Servicios de Tránsito\Documentos\ElectivaUn2\resultados_informe2 

9 Referencias

• Fick, S. E., & Hijmans, R. J. (2017). WorldClim 2: new 1-km spatial resolution climate surfaces for global land areas. International Journal of Climatology, 37(12), 4302–4315. https://doi.org/10.1002/joc.5086
• WorldClim. Historical climate data, version 2.1, periodo 1970–2000. https://www.worldclim.org/data/worldclim21.html
• Food and Agriculture Organization of the United Nations (FAO). Crop requirements and limitations. https://www.fao.org/4/x5648e/x5648e0e.htm

sessionInfo()

R version 4.5.0 (2025-04-11 ucrt)
Platform: x86_64-w64-mingw32/x64
Running under: Windows 11 x64 (build 26100)

Matrix products: default
  LAPACK version 3.12.1

locale:
[1] LC_COLLATE=Spanish_Colombia.utf8  LC_CTYPE=Spanish_Colombia.utf8    LC_MONETARY=Spanish_Colombia.utf8
[4] LC_NUMERIC=C                      LC_TIME=Spanish_Colombia.utf8    

time zone: America/Bogota
tzcode source: internal

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] htmltools_0.5.9   patchwork_1.3.2   kableExtra_1.4.1  knitr_1.52        scales_1.4.0     
 [6] viridis_0.6.5     viridisLite_0.4.3 ggplot2_4.0.3     tidyr_1.3.2       dplyr_1.2.0      
[11] tmaptools_3.3     leaflet_2.2.3     raster_3.6-32     sp_2.2-0          sf_1.1-3         
[16] geodata_0.6-9     terra_1.9-50     

loaded via a namespace (and not attached):
 [1] gtable_0.3.6       xfun_0.56          bslib_0.10.0       htmlwidgets_1.6.4  lattice_0.22-6    
 [6] vctrs_0.7.1        tools_4.5.0        crosstalk_1.2.2    generics_0.1.4     parallel_4.5.0    
[11] tibble_3.3.1       proxy_0.4-29       pkgconfig_2.0.3    KernSmooth_2.23-26 RColorBrewer_1.1-3
[16] S7_0.2.1           lifecycle_1.0.5    compiler_4.5.0     farver_2.1.2       stringr_1.6.0     
[21] textshaping_1.0.4  codetools_0.2-20   stars_0.7-3        class_7.3-23       sass_0.4.10       
[26] yaml_2.3.12        pillar_1.11.1      jquerylib_0.1.4    classInt_0.4-11    cachem_1.1.0      
[31] lwgeom_0.2-17      abind_1.4-8        tidyselect_1.2.1   digest_0.6.39      stringi_1.8.7     
[36] purrr_1.2.1        labeling_0.4.3     fastmap_1.2.0      grid_4.5.0         cli_3.6.5         
[41] magrittr_2.0.4     XML_3.99-0.25      e1071_1.7-17       withr_3.0.2        rappdirs_0.3.4    
[46] rmarkdown_2.32     otel_0.2.0         gridExtra_2.3      evaluate_1.0.5     rlang_1.1.7       
[51] Rcpp_1.1.1         glue_1.8.0         DBI_1.3.0          xml2_1.5.2         svglite_2.2.2     
[56] rstudioapi_0.18.0  jsonlite_2.0.0     R6_2.6.1           systemfonts_1.3.1  units_1.0-1