CATIE - Interpretación de puntos de muestreo de manglares

Resumen

FIESTA es una biblioteca de software para análisis de datos de inventarios forestales. En este documento se explica el uso de FIESTA mediante un caso de ejemplo basado en un conjunto de datos de manglares interpretados por el CATIE.

Introducción

FIESTA (Forest Inventory Estimation and Analysis) es una biblioteca de software, desarrollada en el lenguaje de programación R, para el análisis de datos de inventarios forestales basados en muestras. Fue desarrollada por el Programa de Inventario y Análisis Forestal (FIA) del Servicio Forestal del Departamento de Agricultura (USDA) de los Estados Unidos.

En este documento se detalla el análisis, mediante FIESTA, de un conjunto de datos de 43389 puntos de muestreo de manglares de Costa Rica. Estos datos fueron interpretados por el Centro Agronómico Tropical de Investigación y Enseñanza (CATIE).

Para cada punto se interpretaron las variables:

• Uso del suelo.

La interpretación se realizó en varios tiempos, de los cuales se están considerando dos para este análisis:

• t1: 1998.
• t2: 2021.

El documento fue desarrollado en el sistema de publicación técnica y científica Quarto, el cual combina código en R y sus salidas (tablas, gráficos, mapas) con texto en Markdown.

Trabajo previo al análisis

Instalación de software

Sistema base de R y herramientas de desarrollo

Para trabajar con FIESTA, se debe instalar:

• El sistema base de R. R es un lenguaje de programación enfocado en análisis estadístico y visualización de datos.
• La interfaz de desarrollo integrada RStudio Desktop, la cual proporciona un editor de texto y otras herramientas para escribir programas en R y visualizar sus resultados, entre otras facilidades. También pueden emplearse otros ambientes de desarrollo que ofrecen características similares (ej. Visual Studio Code).

Versión mínima de R

De acuerdo con su documentación, la versión mínima de R que requiere FIESTA, a la fecha de escritura de este documento (2023-12-02), es la 4.2.0. Puede consultar la versión de su instalación al ejecutar el siguiente comando en la consola de R:

# Versión de R
R.version.string

[1] "R version 4.3.2 (2023-10-31 ucrt)"

La salida del comando anterior debe indicar que la versión de su instalación de R es mayor o igual a 4.2.0.

Paquetes de R

FIESTA

El paquete FIESTA está disponible en CRAN (Comprehensive R Archive Network), un repositorio en línea que alberga una amplia colección de paquetes y extensiones para el lenguaje de programación R, lo que facilita su instalación y actualización a nuevas versiones.

Para usar FIESTA, debe instalarse primero. Puede utilizar la función install.packages().

# Instalación del paquete FIESTA
install.packages("FIESTA")

El resultado de la instalación puede verificarse al cargar el paquete con la función library().

# Carga del paquete FIESTA
library(FIESTA)

Si el comando anterior no genera ningún mensaje de error, FIESTA debe haberse instalado adecuadamente.

Otros

Además de FIESTA, se recomienda instalar los siguientes paquetes para procesamiento y visualización de datos.

# Paquete para el desarrollo de documentos computacionales
install.packages("rmarkdown")

# Colección de paquetes para análisis de datos
install.packages("tidyverse")

# Estilos para gráficos de tidyverse
install.packages("ggthemes")

# Paquete para limpieza de datos
install.packages("janitor")

# Paquete para tablas interactivas
install.packages("DT")

# Paquete para graficación interactiva
install.packages("plotly")

# Paquete para mapas interactivos
install.packages("leaflet")

# Funciones adicionales para leaflet
install.packages("leaflet.extras")

# Funciones adicionales para leaflet
install.packages("leafem")

Luego de instalarlos, debe cargar los paquetes con la función library().

# Carga de paquetes adicionales
library(rmarkdown)
library(tidyverse)
library(ggthemes)
library(janitor)
library(DT)
library(plotly)
library(leaflet)
library(leaflet.extras)
library(leafem)
library(sf) # se instala con FIESTA

Obtención de este repositorio

Este documento, junto con otros similares, forma parte de un repositorio en GitHub, una plataforma en línea para compartir código fuente de aplicaciones, basada en el sistema de control de versiones Git. El repositorio contiene el código fuente del documento y los datos que se utilizan en los ejemplos. Su dirección es https://github.com/mesa-monitoreo-puntos/fiesta.

Puede descargar el repositorio, como un archivo ZIP, de https://github.com/mesa-monitoreo-puntos/fiesta/archive/refs/heads/main.zip

También puede “clonar” el repositorio mediante el comando clone de Git:

# Clonación de este repositorio
git clone https://github.com/mesa-monitoreo-puntos/fiesta.git

Una vez que el repositorio haya sido descargado o clonado, puede abrirse con RStudio o con otra herramienta de desarrollo.

Variables generales

En esta sección se definen algunas variables generales del proceso.

# Ruta a los datos de puntos de muestreo
ARCHIVO_PUNTOS <- "../../datos/catie/ceo-Mangle-1998-2005-2007-2016-2021-sample-data-2022-03-23.csv"

# Colores
COLOR_BOSQUE              = "#266900" # Corine CR - Bosque denso
COLOR_HUMEDAL             = "#00ccf2" # Corine EU - Water courses
COLOR_NO_DETERMINABLE     = "#000000" # Negro
COLOR_PASTOS              = "#ffffa6" # Corine CR - Pastos limpios
COLOR_AGRICULTURA         = "#Becd05" # Corine CR - Mosaico de cultivos
COLOR_OTRAS_CLASES        = "#005acf" # Corine CR - Canales
COLOR_PLANTACION_FORESTAL = "#55ff00" # Corine CR - Plantación forestal

COLOR_NULO                = "#808080" # Gris

# Paleta de colores de tipos de uso de la tierra
COLORES_USO <- 
    c(
        "Bosque"                 = COLOR_BOSQUE,
        "T1-Bosque"              = COLOR_BOSQUE,
        "T2-Bosque"              = COLOR_BOSQUE,
        "Humedal"                = COLOR_HUMEDAL,
        "T1-Humedal"             = COLOR_HUMEDAL,
        "T2-Humedal"             = COLOR_HUMEDAL,   
        "No determinable"        = COLOR_NO_DETERMINABLE,
        "T1-No determinable"     = COLOR_NO_DETERMINABLE,
        "T2-No determinable"     = COLOR_NO_DETERMINABLE,
        "Pastos"                 = COLOR_PASTOS,
        "T1-Pastos"              = COLOR_PASTOS,
        "T1-Pastos"              = COLOR_PASTOS,        
        "Agricultura"            = COLOR_AGRICULTURA,
        "T1-Agricultura"         = COLOR_AGRICULTURA,
        "T2-Agricultura"         = COLOR_AGRICULTURA,
        "Otras clases"           = COLOR_OTRAS_CLASES,
        "T1-Otras clases"        = COLOR_OTRAS_CLASES,
        "T2-Otras clases"        = COLOR_OTRAS_CLASES,
        "Plantación Forestal"    = COLOR_PLANTACION_FORESTAL,
        "T1-Plantación Forestal" = COLOR_PLANTACION_FORESTAL,
        "T2-Plantación Forestal" = COLOR_PLANTACION_FORESTAL,
        "Nulo"                   = COLOR_NULO,
        "T1-Nulo"                = COLOR_NULO,
        "T2-Nulo"                = COLOR_NULO
    )

# Área de estudio en hectáreas

Carga y limpieza de datos

Los datos de los puntos de muestreo se proporcionaron en un archivo CSV. En el siguiente bloque de código se cargan en el dataframe puntos. Los nombres de las columnas se “limpian” para evitar la presencia de números al inicio y otros problemas que dificultan su manejo. Los valores nulos (NA) se convierten a hileras de texto.

# Carga de datos de puntos de muestreo
puntos <- read_delim(ARCHIVO_PUNTOS)

# Cambio de tipo de datos
puntos <- 
    puntos |>
    mutate(
        plotid = as.numeric(plotid),
        sampleid = as.numeric(sampleid)
    )

# Adición de una columna con números consecutivos, como identificador único
puntos <- 
    puntos |>
  mutate(
    id = row_number()
  ) |>
  select(id, everything())

# Limpieza de los nombres de columnas
puntos <- clean_names(puntos)

# Cambio de nombre de columnas
puntos <- puntos |>
  rename(
    plot_id   = plotid,
    sample_id = sampleid,
    t1_uso    = clase_de_uso_1998,
    t2_uso    = clase_de_uso_21
  )

# Eliminación de espacios en blanco al inicio y al final
puntos <- puntos |>
  mutate(
    t1_uso = trimws(t1_uso),
    t2_uso = trimws(t2_uso)
  )

# Reemplazo de valores NA por la hilera "Nulo"
puntos <- 
    puntos |> 
    mutate(
        t1_uso = replace_na(t1_uso, "Nulo"),
        t2_uso = replace_na(t2_uso, "Nulo")
    )

# Coversión de t1_uso y t2_uso a factores
puntos <-
    puntos |>
    mutate(
        t1_uso = factor(t1_uso, levels = unique(t1_uso)),
        t2_uso = factor(t2_uso, levels = unique(t2_uso))        
    )

El dataframe puntos contiene más de 40 columnas, provenientes del archivo CSV. Algunas de las más importantes para efectos de este análisis son:

• plot_id: identificador de la parcela.
• sample_id: identificador del punto de muestreo.
• t1_uso: interpretación de la variable de uso del suelo en t1.
• t2_uso: interpretación de la variable de uso del suelo en t2.

La siguiente tabla muestra las columnas mencionadas.

# Despliegue de los datos de puntos de muestreo en una tabla
puntos |>
    select(plot_id, sample_id, t1_uso, t2_uso) |>
    datatable(
        caption = "Puntos de muestreo de uso de la tierra",
        rownames = FALSE,
        colnames = c(
            "plot_id", "sample_id", 
            "t1_uso", "t2_uso"
        ),
        options = list(
            pageLength = 5,
            language = list(url = '//cdn.datatables.net/plug-ins/1.10.11/i18n/Spanish.json')
        )
    )

Análisis con paquetes básicos de R

En esta sección, se analizan los datos con paquetes de R como los de la colección de Tidyverse, Plotly R para gráficos estadísticos interactivos y Leaflet R para mapas interactivos.

Uso del suelo

Mapa

El siguiente mapa interactivo, implementado con el paquete Leaflet, muestra la distribución espacial de la interpretación de la variable de uso del suelo en t1 y t2.

# Conversión de datos de puntos de muestreo a objeto sf (vectorial de puntos)
geo_puntos <-
    puntos |>
    select(plot_id, sample_id, lon, lat, t1_uso, t2_uso) |>
    st_as_sf(
        coords = c("lon", "lat"),
        crs = 4326
  )

colores <- 
    colorFactor(
        palette = c(
            COLOR_BOSQUE, COLOR_HUMEDAL, COLOR_NO_DETERMINABLE, COLOR_PASTOS, 
      COLOR_AGRICULTURA, COLOR_OTRAS_CLASES, COLOR_PLANTACION_FORESTAL
        ), 
    levels = c(
        "Bosque", "Humedal", "No determinable", "Pastos",
        "Agricultura", "Otras clases", "Plantación Forestal"
    )
)

leaflet() |>
    addTiles(group = "OSM") |>
    addProviderTiles(
        provider = providers$Esri.WorldImagery, 
        group = "ESRI World Imagery"
    ) |>
    addProviderTiles(
        provider = providers$CartoDB.DarkMatter,
        group = "Dark Matter"
    ) |>        
    addCircleMarkers(
        data = geo_puntos,
        radius = 2,
        fillColor = ~colores(geo_puntos$t1_uso),
        color = ~colores(geo_puntos$t1_uso),        
        # clusterOptions = markerClusterOptions(),
    popup = paste(
      paste0("<strong>Parcela: </strong>", geo_puntos$plot_id),
      paste0("<strong>Muestra: </strong>", geo_puntos$sample_id),
      paste0("<strong>Uso en t1: </strong>", geo_puntos$t1_uso),
      paste0("<strong>Uso en t2: </strong>", geo_puntos$t2_uso),
      sep = '<br/>'
    ),      
        group = "Uso en t1"
    ) |>
  addLegend(
    position = "bottomleft",    
    pal = colores,
    values = geo_puntos$t1_uso,
    title = "Uso en t1",
    group = "Uso en t1"    
  ) |>      
    addCircleMarkers(
        data = geo_puntos,
        radius = 2,
        fillColor = ~colores(geo_puntos$t2_uso),
        color = ~colores(geo_puntos$t2_uso),        
        # clusterOptions = markerClusterOptions(),
    popup = paste(
      paste0("<strong>Parcela: </strong>", geo_puntos$plot_id),
      paste0("<strong>Muestra: </strong>", geo_puntos$sample_id),
      paste0("<strong>Uso en t1: </strong>", geo_puntos$t1_uso),
      paste0("<strong>Uso en t2: </strong>", geo_puntos$t2_uso),
      sep = '<br/>'
    ),      
        group = "Uso en t2"
    ) |>    
  addLegend(
    position = "bottomleft",    
    pal = colores,
    values = geo_puntos$t2_uso,
    title = "Uso en t2",
    group = "Uso en t2"    
  ) |>  
    addLayersControl(
        baseGroups = c("OSM", "ESRI World Imagery", "Dark Matter"),
        overlayGroups = c("Uso en t1", "Uso en t2"),
    )

Distribución de puntos de muestreo en clases de uso

La distribución de los puntos de muestreo en clases de uso del suelo, en t1 y t2, se muestra en un gráfico de barras y en un gráfico de pastel, ambos elaborados con el paquete Plotly.

t1

# Gráfico de barras ggplot2
grafico_barras <-
    puntos |>
    ggplot(aes(x = fct_infreq(t1_uso), fill = fct_infreq(t1_uso))) +
    geom_bar(
    aes(
      text = paste0(
        "Cantidad de puntos: ", after_stat(count)
      )     
    )
    ) +
    scale_fill_manual(values = COLORES_USO) +
    xlab("Uso") +
    ylab("Cantidad de puntos") +
    labs(fill = "Uso") +
    theme_clean() +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    legend.title = element_blank(),
    legend.box.background = element_blank(),
    legend.key = element_blank()
  )

# Gráfico de barras plotly
ggplotly(grafico_barras, tooltip = "text") |> 
  config(locale = 'es') |>
  layout(
    showlegend = TRUE,
    legend = list(
      title = list(text = ''),
      bordercolor = 'rgba(0,0,0,0)' # Establece el color del borde a transparente
    )
  ) 

# Total de puntos
conteo <- puntos |>
    count(t1_uso)

# Asegurarse de que los colores se asignen en el orden correcto
colores_para_plotly <- COLORES_USO[as.character(conteo$t1_uso)]

# Gráfico de pastel plotly
plot_ly(
    conteo, 
    labels = ~ t1_uso, 
    values = ~ n,
    marker = list(colors = colores_para_plotly)
) |>
add_pie()

t2

# Gráfico de barras ggplot2
grafico_barras <-
    puntos |>
    ggplot(aes(x = fct_infreq(t2_uso), fill = fct_infreq(t2_uso))) +
    geom_bar(
    aes(
      text = paste0(
        "Cantidad de puntos: ", after_stat(count)
      )     
    )
    ) +
    scale_fill_manual(values = COLORES_USO) +
    xlab("Uso") +
    ylab("Cantidad de puntos") +
    labs(fill = "Uso") +
    theme_clean() +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    legend.title = element_blank(),
    legend.box.background = element_blank(),
    legend.key = element_blank()
  )

# Gráfico de barras plotly
ggplotly(grafico_barras, tooltip = "text") |> 
  config(locale = 'es') |>
  layout(
    showlegend = TRUE,
    legend = list(
      title = list(text = ''),
      bordercolor = 'rgba(0,0,0,0)' # Establece el color del borde a transparente
    )
  ) 

# Total de puntos
conteo <- puntos |>
    count(t2_uso)

# Asegurarse de que los colores se asignen en el orden correcto
colores_para_plotly <- COLORES_USO[as.character(conteo$t2_uso)]

# Gráfico de pastel plotly
plot_ly(
    conteo, 
    labels = ~ t2_uso, 
    values = ~ n,
    marker = list(colors = colores_para_plotly)
) |>
add_pie()

Análisis con el paquete FIESTA

En esta sección se analizan los datos con el paquete FIESTA y su módulo de inventarios basados en fotografías (Photo-Based). Este módulo calcula estimaciones de población y errores de muestreo asociados. A diferencia de los estimadores tradicionales del libro verde de FIA (utilizados en otros módulos de FIESTA), que se construyeron basándose en el paradigma de muestreo finito utilizando parcelas de muestra con área distinta, los estimadores basados en fotos se construyeron en el contexto del paradigma de muestreo infinito, junto con el concepto de una región de soporte. FIESTA incluye estimadores no proporcionales para estimaciones de área y cobertura porcentual por dominio, y estimadores de media de cocientes para estimaciones de área y cobertura porcentual dentro del dominio, y soporta la post-estratificación para reducir la varianza.

Se calcula el porcentaje de cada clase de uso del suelo en cada parcela y el porcentaje total de cada clase en el área de estudio, tanto en t1 como en t2. Por último, se analiza el flujo de cambios de t1 a t2.

Uso del suelo

Porcentaje de usos por parcela

Se utiliza la función datPBpnt2pct para transponer los datos de puntos a porcentajes en cada parcela de monitoreo.

t1

# Transposición de puntos a porcentajes de usos por parcela en t1
porcentaje_uso_parcela_t1 <- 
    datPBpnt2pct(puntos, uniqueid="plot_id", tvar="t1_uso")

# Despliegue de la transposición de puntos a porcentajes de uso por parcela en t1
porcentaje_uso_parcela_t1 |>
  datatable(
    caption = "Porcentaje de usos por parcela en t1",
    rownames = FALSE,
    options = list(
      pageLength = 5,
      language = list(url = '//cdn.datatables.net/plug-ins/1.10.11/i18n/Spanish.json')
    )
  )

t2

# Transposición de puntos a porcentajes de uso por parcela en t2
porcentaje_uso_parcela_t2 <- 
    datPBpnt2pct(puntos, uniqueid="plot_id", tvar="t2_uso")

# Despliegue de la transposición de puntos a porcentajes de usos por parcela en t2
porcentaje_uso_parcela_t2 |>
  datatable(
    caption = "Porcentaje de usos por parcela en t2",
    rownames = FALSE,
    options = list(
      pageLength = 5,
      language = list(url = '//cdn.datatables.net/plug-ins/1.10.11/i18n/Spanish.json')
    )
  )

Porcentaje total de usos

La función modPBpop() genera datos de poblaciones, los cuales se utilizan posteriormente para calcular el porcentaje total (en todas las parcelas) en t1 y t2.

# Generación de datos de poblaciones
PBpopdat <- modPBpop(pnt = puntos,
    pltassgnid = "plot_id",
    pntid = "sample_id")

t1

# Estimación de la distribución de usos en t1
LCt1 <- modPB(PBpopdat = PBpopdat, rowvar = "t1_uso")

results.LCt1 <- LCt1$est

# Despliegue de la distribución de usos en t1
results.LCt1 |>
  datatable(
    caption = "Distribución de usos en t1",
    rownames = FALSE,
    options = list(
      pageLength = 10,
      language = list(url = '//cdn.datatables.net/plug-ins/1.10.11/i18n/Spanish.json')
    )
  )

t2

# Estimación de la distribución de usos en t2
LCt2 <- modPB(PBpopdat = PBpopdat, rowvar = "t2_uso")

results.LCt2 <- LCt2$est

# Despliegue de la distribución de usos en t2
results.LCt2 |>
  datatable(
    caption = "Distribución de usos en t2",
    rownames = FALSE,
    options = list(
      pageLength = 10,
      language = list(url = '//cdn.datatables.net/plug-ins/1.10.11/i18n/Spanish.json')
    )
  )

Cambios porcentuales en usos de t1 a t2

PBpoparea <- 
    modPBpop(
      pnt = puntos,
      pltassgnid = "plot_id",
      pntid = "sample_id"
  )

# Uso en t1 vs uso en t2
usoT1vT2 <- 
    modPB(
        PBpopdat = PBpoparea,
      rowvar = "t1_uso",
      colvar = "t2_uso"
    )

# Despliegue de datos en una tabla
usoT1vT2$est |>
  datatable(
    caption = "Uso en t1 vs uso en t2",
    rownames = FALSE,
    options = list(
      pageLength = 10,
      language = list(url = '//cdn.datatables.net/plug-ins/1.10.11/i18n/Spanish.json')
    )
  )

Parte de la información de la tabla anterior se muestra en el siguiente gráfico de líneas.

LCt1.extract <- data.frame(LCt1$est$t1_uso, LCt1$est$Estimate)

# Eliminación de la fila "Total"
LCt1.extract <- 
    LCt1.extract |>
  filter(LCt1.est.t1_uso != "Total")

Year1998 <- as.data.frame(c(1998, 1998, 1998, 1998, 1998, 1998, 1998))
names(Year1998)[1] <- "Year"
LCt1.extract <- data.frame(LCt1.extract, Year1998)
names(LCt1.extract)[1] <- "Variable"
names(LCt1.extract)[2] <- "Value"


LCt2.extract <- data.frame(LCt2$est$t2_uso, LCt2$est$Estimate)

# Eliminación de la fila "Total"
LCt2.extract <- 
    LCt2.extract |>
  filter(LCt2.est.t2_uso != "Total")

Year2021 <- as.data.frame(c(2021, 2021, 2021, 2021, 2021, 2021, 2021))
names(Year2021)[1] <- "Year"
LCt2.extract <- data.frame(LCt2.extract, Year2021)
names(LCt2.extract)[1] <- "Variable"
names(LCt2.extract)[2] <- "Value"

LCt.1y2 <- rbind(LCt2.extract, LCt1.extract)

LCt.1y2$Value <- as.numeric(LCt.1y2$Value)

# Gráfico de líneas en ggplot2
grafico_lineas <-
    ggplot(data = LCt.1y2, aes(
        x = Year,
        y = Value,
        color = factor(Variable)
    )) +
    geom_line(size = 1) +
    scale_x_continuous(limits = c(1998, 2021), breaks = c(1998, 2021)) +
    scale_y_continuous(limits = c(0, 100), breaks = seq(0, 100, 10)) +
    xlab("Año de análisis") + 
    ylab ("Porcentaje de uso") +
  scale_colour_manual(values = 
    c("Bosque"              = COLOR_BOSQUE,
      "Humedal"             = COLOR_HUMEDAL,
      "No determinable"     = COLOR_NO_DETERMINABLE,
        "Pastos"              = COLOR_PASTOS,
        "Agricultura"         = COLOR_AGRICULTURA,
        "Otras clases"        = COLOR_OTRAS_CLASES,
        "Plantación Forestal" = COLOR_PLANTACION_FORESTAL,
        "Nulo"                = COLOR_NULO
    )
  ) +
    labs(colour = "Uso del suelo") +
    theme_clean() +
    theme(legend.position = "top")

ggplotly(grafico_lineas) |> 
  config(locale = 'es') |>
  layout(
    showlegend = TRUE,
    legend = list(
      title = list(text = ''),
      bordercolor = 'rgba(0,0,0,0)'
    )
  ) 

Transiciones de t1 a t2

El siguiente gráfico de Sankey muestra las cantidades de puntos que han cambiado en las clases entre t1 y t2.

puntos_sankey <-
    puntos |>
    select(sample_id, t1_uso, t2_uso)

# write_csv(puntos_sankey, "puntos_sankey.csv")
# puntos_sankey <- read_csv("puntos_sankey.csv")

# El prefijo es para diferenciar t1 y t2 en el gráfico de Sankey
puntos_sankey <- 
    puntos_sankey |>
  mutate(
    t1_uso = paste0("T1-", t1_uso),
    t2_uso = paste0("T2-", t2_uso)
  )

label <- c("T1-Bosque", "T1-Humedal", "T1-No determinable", "T1-Pastos",
                     "T1-Agricultura", "T1-Otras clases", "T1-Plantación Forestal",
                     "T2-Bosque", "T2-Humedal", "T2-No determinable", "T2-Pastos",
                     "T2-Agricultura", "T2-Otras clases", "T2-Plantación Forestal"
                    )

color <- c(COLOR_BOSQUE, COLOR_HUMEDAL, COLOR_NO_DETERMINABLE, COLOR_PASTOS, 
           COLOR_AGRICULTURA, COLOR_OTRAS_CLASES, COLOR_PLANTACION_FORESTAL,
           COLOR_BOSQUE, COLOR_HUMEDAL, COLOR_NO_DETERMINABLE, COLOR_PASTOS, 
           COLOR_AGRICULTURA, COLOR_OTRAS_CLASES, COLOR_PLANTACION_FORESTAL
                    )

# Agregar columnas de códigos para usar en el gráfico de Sankey
puntos_sankey <- puntos_sankey |>
  mutate(
    t1_codigo = case_when(
      t1_uso == 'T1-Bosque' ~ 0,
      t1_uso == 'T1-Humedal' ~ 1,
      t1_uso == 'T1-No determinable' ~ 2,
      t1_uso == 'T1-Pastos' ~ 3,
      t1_uso == 'T1-Agricultura' ~ 4,
      t1_uso == 'T1-Otras clases' ~ 5,
      t1_uso == 'T1-Plantación Forestal' ~ 6
    ),
    t2_codigo = case_when(
      t2_uso == 'T2-Bosque' ~ 7,
      t2_uso == 'T2-Humedal' ~ 8,
      t2_uso == 'T2-No determinable' ~ 9,
      t2_uso == 'T2-Pastos' ~ 10,
      t2_uso == 'T2-Agricultura' ~ 11,
      t2_uso == 'T2-Otras clases' ~ 12,
      t2_uso == 'T2-Plantación Forestal' ~ 13
    )
  )

# Contar las combinaciones de t1_codigo y t2_codigo
conteo_combinaciones <- 
    puntos_sankey |>
  count(t1_codigo, t2_codigo)

# Crear las listas separadas
source <- conteo_combinaciones$t1_codigo
target <- conteo_combinaciones$t2_codigo
value <- conteo_combinaciones$n

# Gráfico de Sankey
grafico_sankey <- 
    plot_ly(
    type = "sankey",
    orientation = "h",
    node = list(
      label = label,
      color = color,
      pad = 15,
      thickness = 20,
      line = list(
        color = "black",
        width = 0.5
      ),
      font = list(
        size = 12,
        color = "#000000",
        weight = "bold"
      )
    ),
    link = list(
      source = source,
      target = target,
      value =  value
    )
  )

grafico_sankey <- 
    grafico_sankey |> 
    layout(
    font = list(
      size = 12,
      color = "#000000",
      weight = "bold"
    )
  )

grafico_sankey

Recursos de referencia y consulta

• Forest Inventory Estimation and Analysis (FIESTA) (sitio oficial de FIESTA)
• FIESTA - Photo-Based Module (módulo de inventarios basados en fotografías de FIESTA)