Área de estudio

El área comprende las selvas medianas subcaducifolia (SMSC) y subperennifolia (SMSP) de la península de Yucatán, México (Figura 1). La superficie aproximada del área bajo análisis es de 77 000 km², donde 28 000 km² corresponden a SMSC y 49 000 km² a la SMSP (^{Instituto Nacional de Estadística y Geografía [INEGI], 2013}).

En la SMSC, el clima es cálido subhúmedo (Aw₁) con lluvia en verano (mayo-octubre) y una estación seca (noviembre-abril). La temperatura media anual es 26 °C (^{Dupuy-Rada et al., 2012}). En la SMSP, el clima es cálido subhúmedo con precipitación media anual de 950 mm (julio-octubre) y temperatura media anual de 22 °C (^{Aryal, De Jong, Ochoa-Gaona, Esparza-Olguin, & Mendoza-Vega, 2014}). En la mayor parte de la península de Yucatán, el gradiente de precipitación es seco hacia el noroeste (600 mm) y más húmedo hacia el sureste (1 400 mm) (^{Martínez & Galindo, 2002}).

Figura 1 Área de estudio: selvas mediana subcaducifolia y subperennifolia -serie V, INEGI, 2013- de la península de Yucatán, México. Las franjas rojas indican la ubicación de datos de biomasa estimados a partir de datos LiDAR (^{Ortiz-Reyes et al., 2019}). 

Datos de biomasa aérea

Se emplearon las estimaciones de biomasa aérea reportadas por ^{Ortiz-Reyes et al. (2019}), correspondientes a archivos ráster. Estos autores emplearon un enfoque basado en el área para la estimación de la biomasa mediante la vinculación de mediciones en campo (datos del Inventario Nacional Forestal y de Suelos de México, remedición 2009-2014) con métricas LiDAR, a través del algoritmo random Forest. Las estimaciones de dicho protocolo, frecuentemente, mantienen e incluso superan los estándares de precisión operacional que las obtenidas mediante el inventario tradicional, con un nivel de sesgo aceptable (^{White et al., 2013}). Cada pixel (20 m x 20 m) que contenía un valor de biomasa aérea estimada sobre las franjas se consideraron como “parcelas LiDAR”. Estas aumentan la distribución y el tamaño de muestra de los datos locales y son similares a las estimaciones de las parcelas en campo (^{Wulder et al., 2012}). Posteriormente, los pixeles de las franjas se remuestrearon a una resolución espacial de 30 m empleando el método del vecino más cercano, para hacerlos coincidentes con los pixeles de las imágenes Landsat (^{Cracknell, 1998}).

Datos de imágenes Landsat y climáticos

Seis imágenes Landsat 8 obtenidas por el sensor OLI (Operational Land Imager) se procesaron a nivel de reflectancia en superficie (^{Vermote, Justice, Clavarie, & Franch, 2016}). Las imágenes se descargaron en la base de datos del Servicio Geológico de Estados Unidos (^{USGS, 2017}).

El periodo de búsqueda en las imágenes fue un año (1 de abril del 2013 hasta el 30 de abril del 2014) para establecer cercanía entre las condiciones de las estimaciones previas de biomasa aérea de los datos LiDAR con las imágenes Landsat. Se descargaron aquellas imágenes con nubosidad menor del 21 %, correspondientes a la época de invierno (Cuadro 1). Una máscara de nubes se aplicó a cada escena empleando el filtro Pixel QA (Quality Assessment) (^{Vermote et al., 2016}). El procedimiento se realizó en el programa QGIS versión 3.6 Noosa (^{QGIS, 2019}), empleando el complemento CloudMasking. Las áreas de nubes y sombras se excluyeron de los análisis subsiguientes.

Posteriormente, los índices espectrales NDVI, MSAVI (índice de vegetación ajustado al suelo modificado), SAVI y EVI (índice de vegetación mejorado) se crearon con las imágenes preprocesadas. Los índices de brillo, verdor y humedad se calcularon mediante la transformación Tasseled Cap para aprovechar la información de más bandas, empleando los coeficientes para productos Landsat con reflectancia en superficie señalados por ^{Crist (1985}). Lo anterior se calculó utilizando el paquete raster en R (^{R Development Core Team, 2013}). También se emplearon las bandas espectrales por sí solas (2 a 7) como variables independientes en la estimación de biomasa aérea.

La información climática se obtuvo de la base de datos del ^{WorldClim (2017)} que dispone de temperatura mensual promedio, mínima y máxima, y precipitación para el periodo 1970 a 2000. Se descargaron los datos de temperatura mensual promedio y precipitación mensual de enero, febrero, noviembre y diciembre, además de la temperatura media anual (°C) y la precipitación anual (mm), variables biológicamente significativas (^{Fick & Hijmans, 2017}). Todas estas variables tenían resolución de ~1 km² por lo que fueron remuestreadas a 30 m, usando el método del vecino más cercano, para hacerlas coincidentes con las demás variables. Estas variables se eligieron por su relevancia probada en otros estudios de estimación de parámetros forestales (^{Ahmed, Franklin, Wulder, & White, 2015}; ^{Wilkes et al., 2015}). La lista de variables predictoras procesadas se muestra en el Cuadro 2.

Estimación de biomasa aérea mediante el algoritmo random Forest

A partir de los datos de biomasa aérea, estimados previamente por ^{Ortiz-Reyes et al. (2019}) en transectos con datos LiDAR (más de 300 000 pixeles para cada tipo de vegetación), se seleccionó aleatoriamente y sin remplazo una muestra de 5 000 pixeles por tipo de vegetación para ajustar sendos modelos y estimar la biomasa para toda el área de interés. La muestra de 5 000 pixeles representó los valores de la variable dependiente (biomasa aérea). Los datos climáticos y los datos espectrales registrados en las imágenes Landsat, correspondientes a los mismos pixeles de la muestra aleatoria, representaron a las variables independientes.

El algoritmo random Forest de R (^{R Development Core Team, 2013}) construye un conjunto de árboles de decisión a partir de datos de entrenamiento, mismos que son validados internamente para generar una predicción de la variable de respuesta dados los predictores (^{Cutler, Cutler, & Stevens, 2012}). El algoritmo es de aplicación fácil y con capacidad de procesamiento de grandes bases de datos de manera eficiente, lo que permite ser una opción en estudios regionales (^{Asner & Mascaro, 2014}). Las variables predictoras definitivas se seleccionaron considerando la influencia que cada una representó sobre el error cuadrático medio (MSE) del modelo ajustado. De manera iterativa, el algoritmo evaluó el rendimiento del modelo para cada tipo de vegetación con base en la RMSE (raíz del error cuadrático medio), el número de variables predictoras y el porcentaje de varianza explicado.

Mapeo de biomasa aérea en los dos tipos de selva mediana

Los mapas de biomasa aérea se crearon en el paquete raster de R (^{R Development Core Team, 2013}). Los mapas se produjeron con el modelo generado previamente utilizando los archivos ráster correspondientes a las variables espectrales y climáticas elegidas por el modelo como relevantes para predecir la biomasa aérea.

Cuantificación de incertidumbre

La incertidumbre refiere el nivel de desconocimiento del valor verdadero de un parámetro o de una variable de interés debido a factores múltiples y se puede cuantificar con estimadores estadísticos comúnes como la desviación estándar, el coeficiente de variación (CV) o mediante un intervalo con un nivel de confianza prestablecido (^{Global Observation of Forest and Land Cover Dynamics [GOFC-GOLD]}). Este estudio evaluó y mapeó la incertidumbre de las predicciones de biomasa aérea a través del CV asociado a las estimaciones generadas a nivel de pixel. El cálculo se realizó con el paquete ModelMap de R (^{Freeman, Frescino, & Moisen, 2018}).

Modelos para estimar biomasa aérea en dos tipos de selva mediana

Se ajustó un modelo por separado para la estimación de la biomasa aérea en cada tipo de selva mediante el algoritmo random Forest. Los parámetros principales se indican en el Cuadro 3.

La capacidad predictiva de los modelos está dentro del intervalo reportado en otros estudios para bosques tropicales (R² = 0.50-0.92). En dichos estudios se emplean datos y enfoques similares a los de esta investigación; por ejemplo, ^{Lu et al. (2012}) estimaron la biomasa aérea en la cuenca amazónica mediante un modelo de regresión múltiple y diferenciaron bosques maduros (R² = 0.50) y en sucesión secundaria (R² = 0.76). Los autores señalan que la estimación de biomasa aérea usando imágenes Landsat es dependiente del sitio, debido a la variación en la fenología, el tipo y estructura de la vegetación. En contraste, ^{Basuki, Skidmore, Hussin, y Van Duren (2013}) usaron imágenes obtenidas por un radar de apertura sintética (SAR, por sus siglas en inglés) e imágenes Landsat ETM+ para la estimación de la biomasa aérea en bosques tropicales bajo manejo en Indonesia. A través de modelos de regresión, dichos autores explicaron 75 % de la varianza (RMSE = 78.9 Mg·ha^-1), mientras que en un bosque tropical de Malasia, ^{Phua et al. (2017}) adjudicaron 63 % de la varianza a las métricas LiDAR y 18 % a las variables de textura de la banda verde de Landsat 8 (RMSE = 112.15 Mg·ha^-1). Por su parte, ^{Ghosh y Behera (2018}) estimaron la biomasa aérea de dos especies plantadas en un bosque tropical en India con datos SAR e imágenes Sentinel-2A; la varianza explicada fue 60 % y 71 % (RMSE = 79.45 Mg·ha^-1; 105.02 Mg·ha^-1) con técnicas de aprendizaje autónomo random Forest y Gradient Boosting, respectivamente. A escala regional, ^{Asner y Mascaro (2014}) estimaron la densidad de carbono aéreo en 14 ecorregiones tropicales de cinco países y, mediante el ajuste de modelos de probabilidad máxima no lineal, explicaron 92.3 % de la varianza (RMSE = 17.12 Mg C·ha^-1).

En los casos anteriores, el rendimiento de los métodos de predicción fue superior al obtenido en el presente trabajo; sin embargo, los errores obtenidos también fueron altos (RMSE entre 78.9 Mg·ha^-1 y 112.15 Mg·ha^-1) en comparación con los de este estudio (RMSE = 34.1 Mg·ha^-1 y 26.2 Mg·ha^-1 para la SMSP y SMSC, respectivamente). No obstante, si los resultados se comparan con el estudio regional de ^{Asner y Mascaro (2014}), el error reportado es similar en términos de biomasa aérea.

Con relación a las variables seleccionadas por el algoritmo random Forest (Figura 2), los resultados mostraron que la banda verde (B3fc) y la precipitación anual (Var_Bio12) fueron las más importantes en los modelos para la estimación de biomasa aérea en la SMSP y SMSC, respectivamente. Las variables climáticas prevalecieron para predecir la biomasa aérea en ambos tipos de selvas medianas, apuntando a una influencia de los procesos climáticos sobre la biomasa aérea estimada. Dicha influencia ha sido reconocida en varios trabajos de estimación de biomasa en vegetación tropical (^{Álvarez-Dávila et al., 2017}; ^{Vieilledent et al., 2016}).

Figura 2 Importancia relativa de las variables seleccionadas por el modelo random Forest en la selva mediana subperennifolia (12 variables) y la selva mediana subcaducifolia (15 variables), para la estimación de biomasa aérea. B3fc: reflectancia de la Banda 3 - verde del sensor Landsat 8 OLI; Prec_M02: precipitación del mes 2 (febrero) (mm); bNDVIfc: Índice de vegetación de diferencia normalizada; Var_Bio12: precipitación anual (mm); Prec_M01: precipitación del mes 1 (enero, mm); Var_Bio1: temperatura media anual (°C); TemAv_M2: temperatura promedio del mes 02 (febrero, °C); B2fc: reflectancia de la Banda 2 - azul del sensor Landsat 8 OLI; WetVal: humedad en la transformación Tasseled Cap; TemAv_M11: temperatura promedio del mes 11 (noviembre, °C); B6fc: reflectancia de la Banda 6 - infrarrojo de onda corta 1 del sensor Landsat 8 OLI; B4fc: reflectancia de la Banda 4 - roja del sensor Landsat 8 OLI; Prec_M12: precipitación del mes 12 (diciembre, mm); GreenVal: verdor en la transformación Tasseled Cap; B5fc: reflectancia de la Banda 5 - infrarrojo cercano del sensor Landsat 8 OLI; brighVal: brillo en la transformación Tasseled Cap; B7fc: reflectancia de la Banda 7 - infrarrojo de onda corta 2 del sensor Landsat 8 OLI; TeamAv_12: temperatura promedio del mes 12 (diciembre, °C). % IncMSE: incremento porcentual del error cuadrático medio. 

Respecto a la precipitación, la principal limitante de las selvas secas es el agua en el suelo, lo cual podría sugerir la relevancia de la precipitación anual (Var_Bio12) en el modelo de la SMSC, mientras que los promedios mensuales de precipitación (Prec_M01, Prec_M02) permanecieron en la modelación de ambos tipos de vegetación. En este sentido, ^{Cao et al. (2015}) mencionan que el crecimiento de este tipo de bosques no solo varía con la edad, el tipo de suelo o el antecedente del uso de la tierra, sino también con la precipitación.

El hecho de que la precipitación y la temperatura permanecieran como variables relevantes en los modelos generados podría ser un indicio de la relación entre la cantidad de agua disponible y su interacción con la temperatura para influir sobre los procesos de crecimiento en biomasa. En tal caso, tanto la precipitación como la temperatura estarían ejerciendo un control superior sobre la densidad de biomasa aérea en los bosques tropicales, puesto que ambos varían regionalmente y son dependientes de la escala (^{Álvarez-Dávila et al., 2017}; ^{White & Hood, 2004}). ^{Saatchi et al. (2011}) señalan que la variabilidad espacial de la biomasa aérea depende del clima, de procesos de perturbación y recuperación provocados de forma natural y por el hombre, del tipo de suelo y variaciones en la topografía. Por su parte, ^{Martínez y Galindo (2002}) mencionan que la alta variabilidad espacial y temporal de la precipitación, el sustrato geológico y el desarrollo escaso del suelo fueron factores decisivos en la distribución de la vegetación en un área con caracteristicas similares a la evaluada en este estudio.

Para la predicción de biomasa en la SMSC, el algoritmo random Forest seleccionó las mismas variables espectrales o similares que han mostrado buena capacidad predictiva en bosques de condiciones semejantes. ^{Freitas, Mello, y Cruz (2005}) señalan que el NDVI es un buen indicador de la biomasa aérea en bosques secos y caducifolios. De los índices espectrales, el NDVI fue el único que se mantuvo para ambos modelos, el resto de índices se eliminaron ya que no contribuyeron al rendimiento de los modelos. Las bandas azul, verde e infrarroja fueron componentes similares a los reportados en el trabajo de ^{Foody et al. (2001}), quienes señalan la importancia de considerar todas las bandas útiles del sensor y no solo los índices dependientes de la banda roja. Las bandas verde, roja e infrarroja se mantuvieron como variables explicativas en ambos tipos de vegetación, destacando la banda verde (B3fc). Dicha relevancia también fue reportada por Foody et al. (2001) y ^{Phua et al. (2017}).

Mapeo de biomasa aérea en los dos tipos de selva mediana

La Figura 3 muestra la variación espacial de la biomasa aérea en los dos tipos de selva mediana. La biomasa de la SMSP varió entre 4.0 y 185.7 Mg·ha^-1 (media [ȳ] = 85.2; desviación estándar [s] = 23.2), intervalo que se encuentra dentro de los resultados reportados por otros autores para el mismo tipo de vegetación. ^{Aryal et al. (2014}) reportaron valores de 11.72 a 99.56 Mg C·ha^-1 en bosques secundarios de cuatro años (s = 4.92) y en bosques maduros (s = 20.83), lo cual es similar al intervalo reportado en este estudio. Recientemente, ^{Hernández-Stefanoni et al. (2020}) reportaron 127.5 Mg·ha^-1 de biomasa aérea promedio y un CV menor de 40 %.

En el caso de la SMSC, la biomasa aérea osciló entre 11.7 y 117 Mg·ha^-1 (ȳ = 51.1; s = 17.5). Este valor se encuentra dentro del intervalo reportado por ^{Dupuy-Rada et al. (2012}) para bosques tropicales secos de la península de Yucatán (ȳ = 56 Mg·ha^-1). La similitud podría deberse a que los datos de biomasa aérea de los dos estudios provienen de mosaicos de fragmentos de bosque en edades sucesionales y arreglo espacial diferentes. Para este tipo de vegetación, ^{Hernández-Stefanoni et al. (2020}) reportaron 100.4 Mg·ha^-1 de biomasa aérea promedio y ^{Dai et al. (2014}) estimaron 5.0 a 115.0 Mg C∙ha^-1 con ȳ = 56.6 Mg C∙ha^-1.

Figura 3 Distribución espacial de la biomasa aérea promedio (Mg·ha^-1) en las selvas mediana subperennifolia y subcaducifolia de la península de Yucatán, México. 

Por otra parte, existen estimaciones que reportan cantidades mayores de biomasa aérea que la reportada en el presente estudio; e.g., ^{Hernández-Stefanoni et al. (2014}) reportaron valores medios de biomasa de 109.71 Mg·ha^-1 y 376.77 Mg·ha^-1 para la SMSC y SMSP, respectivamente, cuando utilizaron parcelas de muestreo en campo de 1 000 m². Estos mismos autores reportaron valores medios de biomasa de 147.2 y 270.2 Mg·ha^-1 para la SMSC y SMSP, respectivamente, cuando utilizaron parcelas de muestreo en campo de 400 m² en la misma área de estudio. Esto muestra la complejidad de la comparación de resultados entre estudios de propósito similar, pero que utilizan métodos o enfoques de análisis diferentes, particularmente cuando la dimensión de las áreas bajo análisis es desigual y los elementos paisajísticos son contrastantes como resultado de estados sucesionales espacialmente azarosos, a causa de actividades antrópicas y disturbios naturales (^{Aryal et al., 2014}; ^{Dupuy-Rada et al., 2012}).

Incertidumbre espacial de las predicciones de biomasa aérea

La variabilidad espacial de la incertidumbre de las predicciones se expresó como el porcentaje de CV (Figura 4). Para la SMSP, el CV osciló mayoritariamente entre 25 y 75 % y fue de magnitud mayor para la SMSC (75 a más de 100 %), particularmente en áreas con biomasa aérea promedio baja y aledañas a otras donde se eliminó información, debido al proceso de enmascaramiento para excluir las celdas que contenían nubes. Estos valores de CV, en general, son mayores que los reportados por ^{Hernández-Stefanoni et al. (2020}) para los mismos tipos de vegetación (0 a 75 %, pero mayoritariamente por debajo de 60 %); sin embargo, es importante resaltar el contraste en el tamaño de la superficie analizada en ambos estudios. Los autores antes mencionados analizaron 3 600 km² de cada tipo de vegetación, mientras que en este estudio se analizaron 28 000 km² de SMSC y 49 000 km² de SMSP, por tanto, es de esperar que la variabilidad sea mayor.

Figura 4 Distribución espacial de la incertidumbre (% coeficiente de variación) de biomasa aérea en las selvas medianas subperennifolia y subcaducifolia de la península de Yucatán, México. 

Como la mayoría de las escasas investigaciones, este trabajo utilizó el CV para reportar la variación espacial de la incertidumbre en la estimación de biomasa. Lo anterior pone de manifiesto la importancia de la evaluación de la incertidumbre por componente, para identificar cual aporta mayor error sobre las estimaciones. Por ejemplo, dos componentes que probablemente afectaron los resultados de este trabajo es la discordancia temporal entre los datos en campo y los de sensores remotos, y la no planificación espacial a priori de la colecta de datos con sensores remotos. Otro componente responsable de la incertidumbre observada es el modelo ajustado previamente por ^{Ortiz-Reyes et al. (2019}) para estimar la biomasa aérea en franjas, cuyos datos se utilizaron en este estudio como punto de partida para ajustar un modelo de mayor envergadura a traves del algoritmo random Forest. Por lo tanto, el uso de enfoques que corrijan los errores que cada componente agrega es una tarea pendiente que podría mejorar la precisión de las estimaciones de biomasa aérea.

El análisis realizado deja constancia del esfuerzo actual de cuantificación de biomasa y ofrece puntos de comparación en la ruta de mejora de las metodologías de cuantificación de incertidumbre en ecosistemas forestales complejos. Por otro lado, los resultados representan un intento de estandarización de los reportes de la variacion espacial de incertidumbre como parte importante de la estimación de biomasa aérea forestal.