Revista de Ciencias Ambientales (Trop J Environ Sci) (Julio-Diciembre, 2019). e-ISSN: 2215-3896. Vol 53(2) DOI: https://doi.org/10.15359/rca.53-2.4 Open Acces: www.revistas.una.ac.cr/ambientales e-Mail: revista.ambientales@una.cr Vega-Araya M. y Alvarado-Barrantes R.

---

Análisis de las series de tiempo de variables biofísicas para cuatro ecorregiones de Guanacaste, Costa Rica

Time Series Analysis of Biophysical Variables in Four Ecoregions of Guanacaste, Costa Rica

Mauricio Vega-Arayaa y Ricardo Alvarado-Barrantesb

[Recibido: 30 de noviembre 2018, Aceptado: 4 de marzo 2019, Corregido: 31 de marzo 2019, Publicado: 1 de julio 2019]

1.Introducción

La variación anual cíclica (estacionalidad) de los procesos biofísicos, es un elemento que configura y define los ecosistemas a través de la adaptación a periodos sin lluvia que puede variar entre 5 y 6 meses para ecosistemas tropicales. La región de Guanacaste, Costa Rica, posee ecosistemas que son laboratorios únicos para la ecología estacional, donde un alto porcentaje de las especies son caducifolias como respuesta a este déficit de lluvia.

La estacionalidad define las complicadas relaciones que determinan la edad, la tasa de crecimiento, la reproducción y la mortalidad, las cuales varían a lo largo del año. Las restricciones en la absorción de nutrientes, en la pérdida de hojas y en la reproducción, pueden ser altamente estacionales, lo que permite que los ecosistemas se adapten, de tal forma que algunas actividades fisiológicas solo se realizan durante algunas épocas del año (Varpe, 2017).

El monitoreo espacial y temporal de los parámetros biofísicos de los ecosistemas, proporciona información sobre su estado y respuestas a la estacionalidad (Viña y Gitelson, 2005), pero son los métodos de detección remota, los que permiten realizar el monitoreo espacio-temporal de una forma efectiva y con un bajo costo, lo que permite analizar justamente, la estacionalidad y sus variaciones en el tiempo.

El índice de área foliar (LAI, por sus siglas en inglés), se define como la mitad del área foliar total por unidad de área superficial del suelo, proyectada en el área horizontal local (Jonckheere et al., 2004; Myneni et al., 2001). Por otra parte, el fPAR, por sus siglas en inglés, es la fracción de la radiación fotosintéticamente activa, en otras palabras, la radiación que las plantas absorben en el rango del espectro electromagnético de 0.4 a 0.7 μm (Myneni et al., 2001; Fang et al., 2005a). El índice normalizado de vegetación (NDVI, por sus siglas en inglés), se basa en la diferencia en reflectancia de la banda del infrarrojo cercano con la del rojo del visible (Jackson y Huete, 1991).

El LAI es un atributo ecológico importante que está relacionado con la actividad fotosintética de la vegetación y con la cantidad y la estructura del estrato superior de la vegetación en los ecosistemas. Como tal, del LAI se puede relacionar con la productividad del ecosistema, con el clima y con los estudios ecológicos (min WANG et al., 2007). A partir de este, es posible monitorear el cambio climático y su influencia en el futuro potencial de secuestro de CO2 de los ecosistemas. Por lo anterior, es esencial contar con conocimientos más precisos sobre la ecología de la producción teórica de los diferentes biomas mundiales (humedales, bosques, matorrales o pastizales). El LAI se utiliza ampliamente como variable de entrada para el modelado de la superficie terrestre de los procesos de la biosfera, y especialmente para las predicciones de la producción primaria fotosintética (Chen et al., 1999). Por otra parte, el fPAR es un parámetro utilizado en la detección remota y en el modelado de ecosistemas, ya que representa el intercambio de energía, vapor de agua y CO2 (Fang et al., 2005a). Tanto, LAI como el fPAR son indicadores de productividad del ecosistema, por lo tanto, ambos pertenecen a las variables biofísicas que son útiles para el estudio y relación de la vegetación con la variación estacional ya sea por variabilidad climática como por el cambio climático (Hunt et al., 2004).

Las ecorregiones analizadas en la región de Guanacaste, Costa Rica, presentan áreas heterogéneas donde se combinan zonas productivas y ecosistemas naturales, en particular en la ecorregión de bosques secos del Pacífico, donde se encuentran grandes extensiones de pastos para ganado, arroz, caña de azúcar entre otros. El uso de parámetros biofísicos (LAI, fPAR y NDVI) en áreas extensas, que combina áreas productivas con ecosistemas naturales, hace que los resultados analizados en forma de series de tiempo, sirvan para monitorear el cambio climático y la variabilidad por el efecto del fenómeno El Niño.

Dada la estacionalidad climática, los ecosistemas de la provincia de Guanacaste representan una importante fuente de servicios ecosistémicos, en especial el de regulación hidrológica. Sin embargo, las condiciones climáticas, como los patrones fisiográficos, varían significativamente dentro de Guanacaste y es por ese motivo que se presenta el análisis de cuatro ecorregiones.

El estudio de las series de tiempo es un método que permite mejorar la comprensión de los vínculos entre la productividad ecosistémica, el clima y el cambio de uso de la tierra en los ecosistemas de tierras de Guanacaste.

En este estudio se analizan tres de los parámetros biofísicos que brindan información sobre los ecosistemas como son LAI y la fracción absorbida de la radiación fotosintéticamente activa (fPAR) y el índice normalizado de vegetación (NDVI) para cuatro ecorregiones en la provincia de Guanacaste.

1.1 Productividad primaria como indicador de secuestro de carbono

La productividad de los ecosistemas suele utilizarse para evaluar la cantidad de bienes producidos en una región, como por ejemplo, los alimentos, la madera y servicios ecosistémicos, entre otros, que se proporcionan a la sociedad y se constituye en una de las métricas más importantes utilizadas para evaluar el secuestro de carbono en los ecosistemas terrestres. La producción primaria neta (NPP, por sus siglas en inglés), representa la cantidad neta de carbono que se añade a la biomasa de las plantas como ramas y hojas, y en general, se expresa por unidad de espacio y tiempo (Running et al., 2004). Según Bonan (2016), también se refiere al proceso biofísico mediante el cual la vegetación produce energía química útil neta de la fotosíntesis, es decir, las plantas absorben la energía luminosa y producen carbohidratos a partir del CO2 atmosférico. Esta energía se calcula, usualmente, restando la respiración autotrófica (Ra) de la planta, del total de energía absorbida, también conocida como Productividad Primaria Bruta (GPP, por sus siglas en inglés). La GPP es la tasa por la cual los productores, dentro de un ecosistema, captan y almacenan una determinada cantidad de energía química como biomasa, en un lapso de tiempo dado. El balance principal entre la absorción total de carbono durante la fotosíntesis (GPP) y la pérdida total de carbono durante la respiración (Ra) se explica por las siguientes ecuaciones:

(1)

Donde: NPP es la productividad primaria neta, GPP es la productividad primaria bruta y R es la respiración de las plantas (Luyssaert et al., 2007).

Por su parte la respiración autotrófica (Ra) de las plantas está dada por:

(2)

Donde: Rm es la energía requerida para el mantenimiento y Rc es el crecimiento de las plantas.

NPP es la tasa por la cual todas las plantas en un ecosistema producen energía química neta útil. Es igual a la diferencia entre la velocidad con la cual las plantas en un ecosistema producen energía química útil (GPP) y la velocidad con que utilizan una parte de esta durante la respiración.

La respiración autotrófica (Ra) está compuesta por la suma de los procesos energéticos que se requieren para el correcto funcionamiento de las plantas y su crecimiento (Rc) (Bonan, 2016). Las plantas son productoras netas de biomasa, esto involucra procesos aún más complejos con organismos vivos o heterótrofos que llegan a incidir en el balance de la productividad y usualmente se conoce como respiración heterotrófica (RH). El desequilibrio entre NPP y RH también se conoce como Productividad Neta del Ecosistema (NEP, por sus siglas en inglés) (Luyssaert et al., 2007). Sin embargo, su estimación requiere del uso de variables y mediciones más complejas que, en consecuencia, limitan su uso como indicador final.

Por otra parte, Gower et al., (2001), señalan que el intercambio neto de ecosistemas, el intercambio neto de CO2 entre los ecosistemas terrestres y la atmósfera, es el equilibrio entre dos procesos principales: NPP y Ra.

1.2 El LAI, fPAR y NDVI desde los sensores remotos

La teledetección proporciona la plataforma más factible para observaciones espaciales y temporales continuas de los parámetros biofísicos en las escalas de local a global (Asner et al., 2003; Chen, 1996; Blümel y Reimer, 2009).

El fPAR y el LAI están entre las 50 variables climáticas esenciales reconocidas por el Sistema Mundial de Observación del Clima de las Naciones Unidas (Global Climate Observing System (GCOS, por sus siglas en inglés) (Gobron y Verstraete, 2009).

En términos generales LAI define o se relaciona con el área foliar de la copa o de las partes verdes de la vegetación, mientras que fPAR define la cantidad de radiación solar que ingresa y que es absorbida por la cubierta vegetal.

Se ha encontrado una relación exponencial inversa entre LAI y la intercepción de luz, que es linealmente proporcional a la tasa de producción primaria o NPP. Una parte de esta energía fija es utilizada por los productores primarios para la respiración celular y el mantenimiento de los tejidos existentes, es decir, respiración de crecimiento y de mantenimiento, como se describe en Amthor y Baldocchi (2001). La energía fija restante en forma de masa vegetal (es decir, producto de la fotosíntesis) es la NPP, tal y como se describe en la Ecuación 1. Tanto LAI como fPAR se utilizan, además para calcular la fotosíntesis de la superficie, la evapotranspiración y que a su vez se utilizan para calcular la energía terrestre, el almacenamiento y fijación de carbono, los procesos del ciclo del agua y la biogeoquímica de la vegetación (Myneni, 2015).

Gran cantidad de la energía almacenada por las plantas se emplea en el crecimiento y la reproducción de los productores primarios, mientras que algunos productos de la NPP son consumidos por los herbívoros, en el caso de ecosistemas naturales, o bien, se extrae en forma de alimentos o productos para el ser humano. En las diferentes ecorregiones, la masa de carbono almacenado por unidad de superficie al año (gC/m/año) es utilizada como unidad de medida de productividad, al ser las estimaciones de LAI y fPAR una herramienta útil para predecir la producción primaria foto sintética en GPP y NPP.

Para la estimación de LAI y fPAR en el nivel de paisaje, los datos de sensoramiento remoto se han generado mediante el análisis multiespectral de las firmas de reflectancia de la vegetación (Garrigues et al., 2008). En general, se han utilizado modelos empíricos y semiempíricos para obtener el LAI y fPAR de la reflectancia de los objetos en la superficie, lo que incluye también, índices de vegetación (Garrigues et al., 2008). La mayoría de los algoritmos de estimación de LAI y fPAR se basan en modelos de transferencia de la radiación, que simulan la absorción de radiación y la dispersión en la vegetación (Lotsch et al., 2003). Una revisión detallada de los modelos se puede encontrar en Myneni et al., (1995).

En general, los índices de vegetación (IV) están diseñados para maximizar la sensibilidad a las características de la vegetación y minimizar los factores perturbadores como el ruido atmosférico, la influencia del suelo o la geometría de la iluminación de la toma de las escenas. Las relaciones entre LAI y IV se calibran de manera específica para distintos tipos de vegetación, al utilizar ya sea coexistente in si tu LAI y mediciones de reflectancia o simulaciones de modelos de radiación de dosel (Feng et al., 2006; Myneni et al., 1997). Otras aproximaciones aplican modelos de transferencia radiactiva de reflectancia de superficie y sus diferentes características como el suelo, propiedades ópticas de la hoja y la geometría de la iluminación (Eklundh et al., 2001; Propastin, 2009).

Los usos y coberturas de la tierra (LULC, pos sus siglas en inglés), se han utilizado como auxiliares para la elaboración de los productos de LAI y fPAR para MODIS (Lotsch et al., 2003), al mostrar así, la relación de estas variables con la configuración del paisaje.

La evaluación de la productividad primaria, basada en píxeles, se puede utilizar para mejorar las estimaciones de los parámetros mencionados en la escala continental, con el fin de aumentar la precisión de las predicciones de flujo y proporcionar estimaciones de la NPP. El conocimiento sobre cómo la vegetación cambia con el tiempo y el espacio, explica de cierta forma, cómo la radiación solar reflejada es modificada por la vegetación.

La mayoría de las fuentes de datos de teledetección disponibles durante los últimos años, han sido sensores de satélite con un bajo número de bandas espectrales en las longitudes de onda visibles y NIR (por ejemplo, Advanced Very High Resolution Radiometer, AVHRR por sus siglas en inglés, Landsat). Los parámetros relacionados con la vegetación como LAI o fPAR, solo se puede estimar de manera aproximada. El fPAR es la fracción de radiación activa fotosintética (400-700 ηm), absorbida por la vegetación que, además, es una variable clave en la evaluación de su productividad y las estimaciones de rendimiento y productividad. LAI proporciona un aporte valioso al modelo climático e hidrológico. Al considerar que fPAR muestra una relación lineal positiva con el aumento de NDVI, LAI no está relacionada con NDVI, esto último depende del tipo de vegetación (Friedl et al., 1995). Para estimar LAI y fPAR, a partir de datos de detección remota, los tipos estructurales de dosel deben definirse como aquellos que presentan relaciones NDVI, LAI o fPAR diferentes entre sí. Por lo tanto, muchos esquemas de clasificación, basados en métricas ecológicas, botánicas o funcionales, no son necesariamente adecuados para el cálculo de LAI y fPAR. La utilidad de esta relación depende de la sensibilidad de estas variables a las características de las copas (Friedl et al., 1995). Por ejemplo, los mismos autores encontraron relaciones lineales entre IV, específicamente NDVI y fPAR. Los modelos no lineales dieron un mejor ajuste para la relación IV y LAI. En general, la relación entre fPAR y el NDVI es más fuerte que entre LAI y NDVI. In situ medido fPAR y NDVI normalmente muestran una fuerte relación lineal, lo que sugiere que la covarianza entre fPAR y NDVI es insensible a las variaciones en la distribución del ángulo foliar (LAD Leaf Angle Distribution, por sus siglas en inglés) y heterogeneidad vegetativa (Fensholt et al., 2004). Autores como Iwata et al., (2013) han encontrado que MODIS sobreestima los valores de fPAR, esto podría explicar por qué, aunque existe una correlación lineal entre MODIS, fPAR y NDVI existe una relación lineal, los coeficientes de regresión son más débiles que la relación in situ. Las relaciones fPAR y NDVI no se aplican en la escala mundial (Friedl et al., 1995), por lo que los análisis locales son siempre necesarios.

Todos los algoritmos para obtener estas estimaciones sobre LAI o fPAR son empíricos y se basan en una relación común entre la señal registrada y el respectivo parámetro de las plantas (Fensholt et al., 2004); esta simplificación, a menudo conduce a una baja exactitud en la estimación de estos parámetros y, por último, de la GPP y la NPP.

Las estimaciones de fPAR, basándose solo en este tipo de relación simple, por un lado, no pueden diferenciar entre la absorción por elementos fotosintéticamente activos de las plantas y por otro, la absorción por otras sustancias vegetales y el suelo. Este es el problema principal de los sensores más antiguos como AVHRR; sin embargo, presentan series de tiempo más largas de observaciones en la tierra (Kappas y Propastin, 2012). Los sensores de una nueva generación como MODIS (lanzado en 1999) y MERIS (Medium-Spectral Resolution, Imaging Spectrometer, por sus siglas en inglés), lanzado en 2002, abren nuevas perspectivas. Los de nueva generación y futuros desarrollos de sensores avanzados (por ejemplo, MERIS modificado en SENTINEL-3, EnMap, NPOES), tienen más canales espectrales en las regiones del espectro electromagnético del visible y del infrarrojo cercano en particular. Además, de una mejor resolución espacial, estas plataformas poseen una superior precisión radiométrica. Con el uso y desarrollo de los sensores que puedan dar continuidad a las mediciones, es posible implementar diferentes técnicas para invertir los modelos de la transferencia de radiación en la vegetación, lo que genera series de tiempo. Estos pueden dividirse en modelos de transferencia de radiación de hojas y copas (por ejemplo, PROSPECT y SAIL, como el discutido y propuesto por Jacquemoud et al., (2009) o en un componente de copa simple presentado por (Bicheron, 1999). El modelo de propiedades ópticas de las hojas “PROSPECT” y el de reflectancia bidireccional de las copas “SAIL”, también llamado “PROSAIL” en combinación, se ha utilizado durante aproximadamente dieciocho años para estudiar la reflectancia espectral y direccional del dosel de la vegetación. PROSAIL también se ha utilizado para desarrollar nuevos métodos para la obtención de las propiedades biofísicas de la vegetación, donde se combina la variación espectral de la reflectancia del dosel, que se relaciona, principalmente, con los contenidos bioquímicos de las hojas, con su variación direccional, con la arquitectura del dosel y el contraste suelo-vegetación. Por otra parte, Jacquemoud y Ustin (2003), presentan unas aplicaciones de las estimaciones de las variables biofísicas y estructurales del dosel para aplicaciones en agricultura o ecología, a diferentes escalas.

La imposibilidad de diferenciar entre la absorción por los elementos vegetales fotosintéticamente activos y por la de otras sustancias vegetales y el suelo, es el principal problema de los sensores más antiguos como AVHRR. Los sensores de una nueva generación como MODIS y MERIS de la agencia espacial europea, abren nuevas posibilidades al respecto (Jacquemoud y Ustin, 2003). El modelo PROSAIL también se ha utilizado para desarrollar nuevos métodos para adquirir parámetros de propiedades biofísicas de la vegetación. Este modelo combina la variación espectral de la reflectancia del dosel, que se relaciona especialmente con los contenidos bioquímicos de las hojas, con su variación direccional, y a su vez con la estructura del dosel y la interacción suelo-vegetación (Brown et al., 2000).

Como hemos mencionado, tanto LAI, fPAR y NDVI son factores muy importantes para la caracterización la fenología de las plantas. En general, estas tienen un desarrollo estacional que puede ser cuantificado y analizado como series de tiempo. Dichos análisis son la base de comparación con otros parámetros ambientales como precipitación, temperatura, evapotranspiración, entre otras. Lo anterior es la base para el entendimiento de los procesos ecológicos.

Morisette et al., (2009), presentan esquemáticamente las complejas relaciones entre la fenología y procesos de corto plazo como variaciones del tiempo atmosférico o bien procesos de largo plazo como los cambios climáticos, donde el estudio de la fenología es cada vez más relevante para un mejor entendimiento de los problemas ambientales.

2.Metodología

2.1 La ecorregiones de Guanacaste

El área de estudio se ubica en la provincia de Guanacaste, al noroeste de Costa Rica, América Central. De esta provincia, se extrajeron cuatro ecorregiones de las cinco que contiene. En la Figura 1 se muestran las ecorregiones seleccionadas, que según Olson y Dinerstein (2002) se denominan: a) bosques húmedos estacionales, que se dividen en dos ecorregiones: (bhe1 y bhe2); b) bosques secos del Pacífico (bsp); c) bosques montanos de Talamanca (bmt).

La ecorregión de bosques secos del Pacífico (bsp), se extiende a lo largo de la costa del Pacífico, corresponde a un hábitat tropical que alberga importantes especies de plantas y animales, así como un importante grado de endemismo. Esta se encuentra totalmente fragmentada, representada en menos del 2 % del hábitat original, está amenazada por fuertes presiones provenientes del ser humano, como por ejemplo la ganadería, la quema, la expansión agrícola y las operaciones de caza. Los bosques tropicales secos son ahora mucho más raros que los bosques pluviales tropicales.

En esta ecorregión, la estación seca se prolonga de 5 a 6 meses, con una precipitación promedio anual de entre 1 000 y 2 500 mm y un patrón, generalmente bimodal de lluvia, con un período seco más corto conocido como veranillo. La estación seca se extiende desde finales de noviembre hasta inicios de mayo. Los vientos predominantes son los alisios que, generalmente tienen una dirección de noreste y este hacia el suroeste y sur. La mayor parte de ella se ve influenciada por el sistema montañoso que va del noroeste al sureste en lo que se denomina la cordillera de Guanacaste.

La ecorregión de bosques húmedos estacionales se caracteriza porque muchos de sus árboles son caducifolios, es decir que pierden sus hojas durante la estación seca, esto también es común, a la ecorregión de bosques secos del Pacífico. La flora está más adaptada a la ausencia de lluvia de la estación seca; los animales también se acoplan a esta fluctuación entre húmedo y seco.

Figura 1. Principales ecorregiones en Guanacaste para el análisis de LAI, fPAR y NDVI.

Los bosques montanos de Talamanca (bmt), forman una ecorregión de bosque de montaña que pertenece al bioma de los bosques húmedos latifoliados tropicales y subtropicales, según la definición del Fondo Mundial para la Naturaleza (Programa REDD/CCAD-GIZ-SINAC, 2015). Se encuentra en Costa Rica y Panamá, e incluye la cordillera de Guanacaste, cordillera de Tilarán, cordillera Volcánica Central, y la cordillera de Talamanca (Programa REDD/CCAD-GIZ-SINAC, 2015). En Guanacaste la ecorregión posee un área de 124 472 hectáreas y comprende bosques con una altitud de 800 a 1 900 msnm. Propiamente en la cordillera de Talamanca en Costa Rica, esta ecorregión puede alcanzar elevaciones de más de 3 000 msnm. Para su análisis se usó solamente el área que corresponde a la provincia de Guanacaste, Costa Rica, tal y como se muestra en la Figura 1. Aproximadamente el 75 % de la cobertura forestal original permanece intacta (Programa REDD/CCAD-GIZ-SINAC, 2015), y el 40 % de la ecorregión se encuentra dentro de áreas protegidas, resaltándose los parques nacionales de Rincón de la Vieja, Guanacaste (parte alta), el Volcán Tenorio y la Zona Protectora de Miravalles.

2.2 Análisis efectuados para las cuatro ecorregiones

El espectrorradiómetro de imágenes de resolución moderada (MODIS, por sus siglas en inglés), posee 36 bandas espectrales con diferentes resoluciones espaciales, por ejemplo, dos con 2 bandas que poseen un tamaño de píxel de 250 m, cinco con 5 bandas de 500 m y el resto (29) poseen un píxel de 1 km. El satélite Terra fue lanzado en 1999, el Aqua en 2002, de tal manera que hay un gran número de productos para el estudio de la atmósfera, tierra y agua, los cuales son fácilmente disponibles.

El Departamento de Geografía de la Universidad de Boston, entrega los productos de LAI globales, basados en la Fracción Terra/Aqua-MODIS de productos de radiación activa fotosintética absorbida (fPAR)(MCD15A3H V6). Para una descripción detallada del producto MOD15A2, véase Myneni (2015); Knyazikhin et al., (1998). El último producto MCD15A3H (Colección 6), reemplazó los anteriores LAI y fPAR referenciado de 6 biomas con un nuevo producto referenciado a 8 biomas. Las clases de vegetación con hojas anchas y coníferas fueron divididas en subclases deciduas y de hoja perenne. Los nuevos productos presentan mejoras en el algoritmo de cálculo con lo cual se mejora la calidad de LAI y fPAR. Los valores de fPAR y LAI son sin dimensión y poseen un tamaño de píxel de 500 m. El fPAR se calcula para el rango del espectro electromagnético del visible (400-700 ηm).

Los productos MODIS: LAI, fPAR y NDVI están disponibles en https://urs.earthdata.nasa.gov/(Didan, 2015; Myneni, 2015), así como desde la plataforma Earth Engine. En general, para la generación de las imágenes, los proveedores elaboran los algoritmos de cálculo donde se elige el mejor píxel disponible de todas las adquisiciones de ambos sensores MODIS (de los satélites Aqua y Terra) de la NASA dentro de un período de 4 días.

De los productos MOD13A1 y MYD13A1 de MODIS, se extrae el NDVI, el cual se basa en las propiedades de absorción de la vegetación en las bandas roja e infrarroja, el cual se expresa matemáticamente de la siguiente forma: (Price, 1994; Tucker, 1979):

(3)

Donde: ρROJO es la fracción del espectro electromagnético en el rango de 0.6-0.7 μm, ρINR es la porción del infrarrojo cercano en el rango del espectro electromagnético de 0.75-1.5 μm para las respectivas bandas de MODIS.

Los productos globales de LAI, fPAR y NDVI están disponibles como compuestos mensuales para el período comprendido entre febrero de 2000 y el último mes calendario. En este estudio se usó un tamaño de píxel de 500 m para los tres productos. Se usó la plataforma Google Earth Engine (Gorelick et al., 2017) para la obtención, procesamiento y reducción de las series temporales y para efectuar los procesos que se detallan más adelante.

Se usaron las identidades trigonométricas de seno y coseno para definir los ciclos de las series de tiempo de las tres variables, al considerar la periodicidad, es decir, se definió un ciclo como un periodo completo de una función de seno o bien de coseno en el intervalo de tiempo, como sigue (Shumway y Stoffer, 2017):

(4)

Donde: t=0,±1,±2,…, y ω es un índice de frecuencia (ciclos por unidad de tiempo), A es la amplitud en la función y φ es la fase, que determina el punto de inicio de la función.

Definiendo en la Ecuación 4 la identidad trigonométricac como sigue (Shumway y Stoffer, 2017):

(5)

Donde: se suponen como variables independientes aleatorias de distribución normal. Esto último es media 0 y varianza.

El propósito de este análisis es el de ajustar un modelo de primer orden a las series temporales de datos MODIS. Se utiliza la regresión lineal como reductor de la serie temporal y para estimar términos de seno, coseno y la tendencia de los datos.

Para ajustar este modelo a la serie temporal de LAI, fPAR y NDVI, se definió ω=1 (un ciclo por unidad de tiempo) y posteriormente se usó la técnica de regresión de mínimos cuadrados, de la forma:

(6)

donde: β0 es una constante, β1 es el coeficiente de la variable temporal, β2 =Acos(φ) y β3 =Asin(φ) y A = (β22 + β32 )0,5 y φ = atan (β3 /β2). Por último, e1 es el término de error aleatorio, el cual sigue una distribución normal.

En general, la representación de series de tiempo complejas en la forma de entidades trigonométricas como seno y coseno, se denomina análisis de “Fourier” (Davis, 2003).

Posteriormente a la obtención de las series temporales por ecorregión, se obtuvieron componentes estacionales: tendencia (trend), estacionalidad (seasonal) y componentes irregulares o residuos (remainder), TSR por sus siglas en inglés (Cleveland et al., 1990), al usar el lenguaje de programación R (R Core Team, 2018). Para la correcta descomposición de las series, se graficó cada una de las series y se identificó si cada serie es aditiva o multiplicativa. TSR es un método muy versátil y robusto para descomponer series de tiempo y utiliza el método “LOESS”, que se usa para estimar relaciones no lineales (Cleveland et al., 1990). Los dos grandes parámetros que se configuran al usar TSR son la ventana de tendencias (t.window) y la ventana de estacionalidad (s.window). Para este estudio se configuraron las ventanas en valore de t.window en 10 y s.window en “periodic”. Estos componentes controlan qué tan rápido cambian los componentes de tendencia y estacionalidad.

La tendencia marca la direccional en el largo plazo de la serie de tiempo, en segundo lugar, la estacionalidad marca los movimientos regulares (estacionales, por ejemplo) y sistemáticos, de tal forma que es posible detectar fluctuaciones adicionales o irregulares, no sistemáticas y de corto o mediano plazo (Künzer et al., 2015).

Se realizaron cuatro pruebas para evaluar la estacionalidad de las tres series de tiempo, esto con el fin de escoger el tipo de modelo (aditivo o multiplicativo) de descomposición de las series. La primera, es la denominada prueba Phillips-Perron (Phillips y Perron, 1988), la segunda es la Augmented Dickey-Fuller, por su nombre en inglés (Fuller, 1996). La tercera fue la Kwiatkowski-Philips-Schmidt-Shin (Hobijn et al., 1998). La última prueba fue la Ljung-Box (Box y Pierce, 1970). Las series de tiempo de LAI, fPAR y NDVI se analizaron en sus componentes y según el modelo aditivo, ya que según las cuatro pruebas de estacionalidad dieron que las series se comportan de forma aditiva, ya que en todas las series de LAI, fPAR y NDVI las variaciones estacionales son relativamente constantes, es decir, los valores no incrementan sistemáticamente conforme incrementa el tiempo. Matemáticamente el modelo aditivo de descomposición de las series es: (Brockwell y Davis, 2002):

(7)

donde: Xt es la serie de tiempo del modelo aditivo, ya sea LAI, fPAR y NDVI, mt es la componente de tendencia de la serie, st es el componente estacional y son los remanentes que son atribuido a efectos aleatorios o fluctuaciones aleatorias.

Con la descomposición de las series se busca (Brockwell y Davis, 2002; Forkel et al., 2013; Hyndman y Athanasopoulos, 2013):

•Comprobar si existe una tendencia clara en las series de tiempo (fuerte o débil)

•Identificar la estacionalidad, si es un patrón regular y si hay picos y puntos bajos (baja, media o alta)

•Identificar las variaciones de tiempo corto (baja, media o alta)

•Verificar si hay ciclos de picos o puntos bajos que no corresponden a la estacionalidad

•Identificar los valores atípicos

•Comprobar si existen variaciones en la tendencia (gradual o abrupta) y si son constantes o no constantes en el tiempo

•Ver si existen ciclos, que son diferentes de la estacionalidad (tiempo)

Se graficaron, por sub ecorregión y por variable, la tendencias, estacionalidad y remanentes de las series de LAI, fPAR y NDVI, según Hyndman y Khandakar (2008). Por otra parte, se usó la función “monthplot” (R Core Team, 2018), la cual extrae las subseries mensuales de la serie total. Los resultados se presentan en las Figuras 6, 8 y 10.

3. Resultados

En la Figura 2, se presenta la media de los valores de fPAR para el producto MCD15A3H, en la versión 6 de MODIS al nivel 4 para octubre de 2016 y octubre de 2017. Es claro que la región bmt presenta una evidente diferencia en los valores de fPAR en relación con el resto de las ecorregiones.

Figura 2. Distribución espacial de fPAR (sin dimensión) para las cuatro ecorregiones en Guanacaste.

En la Figura 3 se presenta la media de los valores para LAI para el producto MCD15A3H, en la versión 6 de MODIS en el nivel 4 para octubre de 2016 y octubre de 2017. Los valores de LAI al igual que fPAR, son sin dimensión y poseen un tamaño de píxel de 500m. De la Figura 3, resalta que para la ecorregión bsp existen dos bloques al centro y sur de esta ella, que presentan valores altos de LAI. Estos bloques son parte de dos áreas silvestres protegidas a saber, parque nacional Palo Verde y la reserva biológica Lomas Barbudal.

Figura 3. Distribución de LAI (sin dimensión) para las cuatro ecorregiones en Guanacaste.

En el Apéndice 1, se muestra a manera de ejemplo la imagen de “amplitud” para el NDVI. Tanto la amplitud como la fase (Apéndice 2), se utilizan para descomponer las series de tiempo en el análisis de Fourier.

Figura 4. Distribución de NDVI (sin dimensión) para las cuatro ecorregiones en Guanacaste.

En el Cuadro 1, se presenta un resumen con los valores medios, mínimo y máximo por ecorregión. Estos representan, por ejemplo, para la media, los valores individuales de todos los píxeles de la respectiva ecorregión y en toda la serie de tiempo. Para el caso del NDVI los valores normalmente se reportan entre -1 y 1, sin embargo, para MODIS los productos se escalan en un factor de 10 000.

Cuadro 1. Valores medios mínimos y máximos por ecorregión en Guanacaste para fPAR, LAI y NDVI para la serie temporal de julio de 2002 a octubre de 2017

Variable	bhe 1	bsp	bmt	bhe 2
fPAR media	57.9	54.7	44.6	57.9
fPAR min	29	0	20.1	29
fPAR max	83.6	85.2	70.1	83.6
LAI media	25.9	19.7	22.4	25.9
LAI min	7.1	0	6.4	7.1
LAI max	53.8	56.5	42.4	53.8
NDVI media	7 630	6 989	6 281	7 517
NDVI min	-1 786	-1 693	3 252	-455
NDVI max	9 006	8 875	8 804	9 138

En Apéndice 1, se presenta la distribución media de los valores de amplitud derivados del NDVI.

En la Figura 5, se presenta la serie de tiempo de NDVI descompuesta y por cada una de las cuatro ecorregiones. Esta serie es estacionaria en todas las pruebas. Por otra parte, se puede notar la similitud en la estacionalidad de las ecorregiones bhe1 y bhe2 y contrasta mucho con la estacionalidad de bmt que posee una variación interanual muy marcada. Para el periodo entre el 2007 y mediados de 2010, estas ecorregiones bhe1 y bhe2 tienden a disminuir el NDVI sistemáticamente, para caer en el pico más bajo de la serie en 2011. Por otro lado, el bsp presenta una clara tendencia a la disminución a lo largo de toda la serie con su pico más bajo en 2015.

Figura 5. Serie de tiempo (data) para el periodo 2002-2017, con estacionalidad (seasonal), tendencia (trend) y componentes irregulares (remainder) para el NDVI para las cuatro ecorregiones en Guanacaste.

En las Figuras 6, 8 y 10, se presentan las variaciones mensuales por variable, por ecorregión para las series de tiempo. La línea azul en cada mes representa la media de la serie para ese mes y la línea en negro es la serie completa para ese mes.

Específicamente en la Figura 6, se aprecia que para bhe1 y bhe2 en la variable NDVI, estas dos ecorregiones poseen valores muy similares en cuanto al comportamiento medio de todos los meses, siendo marzo y abril los meses que más varían en la serie. Sin embargo, octubre es un mes que también varía mucho, en contraposición a los meses de junio y julio que varían poco. Resalta también el descenso de la media (línea azul) en la serie de bhe1 y bhe2, a partir de julio y con un mínimo en septiembre para luego incrementar sus valores medios. La ecorregión bsp, presenta una disminución de la media entre los meses solo que a partir de mayo y hasta septiembre, esto coincide también con el inicio de las lluvias (ver Apéndice 5).

La ecorregión bmt es sumamente variable entre los meses de la serie, por ejemplo el mes de mayo en el 2009 contiene el valor menor (3 000) para este mes y en 2013 presenta el máximo (7 000).

Las ecorregiones bhe1, bhe2 y bsp, presentan un descenso en los valores medios de NDVI a partir del mes de julio, lo que podría ser un indicador de saturación (Asner et al., 2003; Chen et al., 2006; Hobbs, 1995; Huete et al., 2002) del NDVI.

Figura 6. NDVI por mes (periodo 2002-2017) para las cuatro ecorregiones en Guanacaste.

Las series de fPAR para las cuatro ecorregiones son estacionarias, según las pruebas realizadas para este fin. La Figura 7 contiene la serie descompuesta para la variable fPAR. Al igual que en la serie de NDVI, la estacionalidad de las ecorregiones bhe1 y bhe2 son muy similares en su variación interanual con un pico anual y un descenso, luego un pequeño ascenso y finalmente, llega al punto más bajo. La ecorregión bmt tiene una estacionalidad con variación interanual muy diferente a bhe1 y bhe2, ya que luego del pico máximo baja hasta el pico mínimo, sube un poco y vuelve a bajar para finalmente, llegar al punto máximo. La estacionalidad de la ecorregión bsp solo posee un leve descenso interanual para pasar del pico mínimo al máximo e iniciar otro ciclo.

Figura 7. Serie de tiempo (periodo 2002-2017) con estacionalidad, tendencia y componentes irregulares para el fPAR para las cuatro ecorregiones en Guanacaste.

El LAI para las cuatro ecorregiones se presenta en la Figura 9. A diferencia de la variable fPAR, las ecorregiones bhe1 y bhe2 difieren un poco en la estacionalidad. La primera diferencia es en magnitud, donde la ecorregión bhe1 muestra más fuerza en lo estacional que bhe2, lo segundo es la variación interanual de bhe1, que es mucho más marcada que bhe2. En lo referente a la tendencia bhe1 y bhe2 son en realidad bastante similares, que tiene un periodo claro de descenso desde 2003 a 2006, luego se hace irregular con un pico mínimo en 2010. Llama la atención un punto de quiebre donde inicia el pico máximo en bhe1 en 2015 y llega al máximo en el tercer trimestre de 2015, luego desciende bruscamente en 2016, donde la anomalía NIÑO 3 (Apéndice 7) presenta los valores positivos (fase cálida) más altos.

La ecorregión bsp presenta una estacionalidad con un leve descenso interanual y claramente varía anualmente con un pico máximo entre octubre y noviembre. El ascenso hacia el pico máximo se da en sincronización con el inicio y descenso de las lluvias. En relación con la tendencia, esta ecorregión presenta dos picos mínimos en 2005 y 2010, sin embargo, es bastante irregular.

Los bosques montanos de Talamanca (bmt) tienen una leve propensión a marcar una estacionalidad interanual, pero muy leve. Resalta que el inicio del descenso desde el pico más alto se da al inicio del cada año para llegar a un punto mínimo hacia la mitad del año. La tendencia del bmt es sumamente irregular que tiene un pico máximo en 2014, el cual resalta en esta curva.

Figura 8. fPAR por mes (periodo 2002-2017) por ecorregión en Guanacaste.

En relación con la variación mensual de fPAR, presentada en la Figura 8, las ecorregiones bhe1, bhe2 presentan un comportamiento similar con un valor medio mínimo en abril y un máximo en julio. Luego de este mes hay un leve descenso. Para estas ecorregiones la variación en la serie para los meses de julio, agosto y septiembre es poca, si la comparamos con los meses de marzo y abril.

La ecorregión bmt presenta una variación de la media mensual muy similar a la variación media mensual de la precipitación, con valores bajos en los meses secos y máximos en agosto, septiembre y octubre. Sin embargo, hay mucha variación en la serie entre los meses de junio a octubre, es decir periodos de dos años, donde los valores de fPAR varían mucho para estos meses.

Los valores medios de la serie de la ecorregión bsp, presenta el pico máximo en junio y el mínimo en octubre. El pico máximo en junio coincide con las ecorregiones bhe1 y bhe2.

En la Figura 9, se presenta la serie para la variable de LAI. La estacionalidad de bhe1 y bhe2, son bastante diferentes, el componente estacional es más fuerte en bhe1, con un periodo de descenso a inicio del año, luego un ascenso que no llega al máximo, para luego bajar nuevamente y al final, ascender hasta el máximo. La ecorregión bhe2 no es tan estacional como bhe1 y el descenso anual no es tan marcado como bhe2.

La ecorregión bsp casi no presenta descenso entre cada año y difiere mucho de la estacionalidad de las otras. En bmt, por ejemplo, hay un leve descenso al inicio del año, para continuar el recorrido hacia el máximo, es decir, el recorrido hacia el máximo se interrumpe, desciende un poco y luego continúa subiendo hacia el máximo. La ecorregión bsp, por su parte, llega casi al máximo y desciende un poco para luego llegar al máximo.

Figura 9. Serie de tiempo (periodo 2002-2017) con estacionalidad, tendencia y componentes irregulares para el LAI para las cuatro ecorregiones en Guanacaste.

Las tendencias de bhe1 y bhe2 son muy similares en su irregularidad. Si se aprecia un periodo entre 2003 y finales de 2005, donde hay un descenso sistemático en el LAI, en las ecorregiones bhe2, bhe2 y bsp. Resalta que el máximo de la ecorregión bhe1, es en el año 2015, justo antes del fuerte evento de Niño de 2015-2016 (ver Apéndice 7), mientras para 2016 desciende.

Figura 10. LAI por mes (periodo 2002-2017) por ecorregión en Guanacaste.

Sobre la variación mensual de las series, se aprecia que las ecorregiones bhe1 y bhe2 son muy similares, los meses secos de diciembre a abril son muy variables. Los meses más lluviosos (agosto-octubre) son constantes en la serie, pero presentan valores bajos en la media mensual de estas dos.

Por otra parte, la ecorregión bmt presenta, en general, mucha variación respecto a la media en todos los meses de las series mensuales. En contraposición con bhe1 y bhe2, los valores medios son en los meses más lluviosos (septiembre y octubre), lo cual se antepone con la ecorregión bsp, donde los valores medios más bajos están en estos dos meses.

En el Apéndice 6 se presentan las líneas de tendencia para precipitación, fPAR, NDVI y LAI para las cuatro ecorregiones, en 2002-2017. En esta figura se puede apreciar, por ejemplo, la similitud de LAI y fPAR en todas las ecorregiones.

4.Discusión

En este artículo se presenta el análisis de variables bioclimáticas (NDVI, LAI y fPAR) en forma de series de tiempo. Esto como primera etapa del estudio de las relaciones entre la precipitación y las variables bioclimáticas, las cuales al ser estimaciones a partir del “verdor” (contenido de clorofila foliar) de las copas o de la parte superior de la vegetación, pueden relacionarse con el área foliar, la cubierta del dosel y la estructura (Huete et al., 2008; Asner et al., 2003) y responden directamente al acceso y disponibilidad de condiciones como agua, luz, temperatura, etc. Por otra parte, y a su vez, el índice de vegetación NDVI se relaciona con variables biofísicas/bioquímicas del dosel, específicas, como LAI y fPAR (Huete et al., 2008). Esta es la base y motivo para analizar la serie de tiempo del NDVI en conjunto con LAI y fPAR y su relación con la precipitación. Otros autores como Choudhury (1987) y Tucker (1979) han relacionado o bien enunciado que LAI y fPAR se relacionan con el contenido de clorofila, fracción verde de la vegetación, con biomasa y los procesos biofísicos del dosel como fotosíntesis y transpiración.

Un aspecto importante para poder analizar estas variables como series de tiempo, es su disponibilidad de procesamiento en la plataforma gratuita de uso público Google Earth Engine. La cual permite extraer para áreas específicas, información sobre diversos productos derivados de diferentes sensores. Un análisis posterior, con el uso de esta plataforma, permitirá integrar las series de AVHRR con las de MODIS. Ya autores como Brown et al., (2006), han sugerido que dada la correlación y similaridad de los sets de datos, es posible reconstruir una serie de tiempo más larga y consistente. Sin embargo, aunque la presentada en este trabajo es de cerca de 15 años, y autores como Baldocchi y Wilson (2001), señalan que, para lograr un monitoreo significativo de la vegetación a largo plazo, los registros deben poseer por lo menos una década, el poder integrar series históricas es sumamente necesario, sobre todo para analizar las series biofísicas con series climáticas.

Este estudio ha utilizado la metodología de descomposición aditiva para el análisis de las series en tendencia y el patrón estacional para las variables. Al utilizar este método y aislar el componente estacional, se facilita el entendimiento de la tendencia de las series de tiempo, al mismo tiempo que se avanza en una exploración que permite identificar la variación en ellas para comparar ecorregiones. El análisis de cuatro ecorregiones para la provincia de Guanacaste ha permitido encontrar similitudes y diferencias. El comparar unidades geográficas definidas (ecorregiones) permite resumir y sintetizar la variabilidad espacial en valores de tendencia central como la media o la mediana, con lo que se puede hacer la comparación no solo entre ecorregiones, sino también monitorear cada ecorregión en el tiempo.

El NDVI como parámetro, que permite interpretar biofísicamente la actividad con la que la vegetación está absorbiendo la fracción del rojo por un lado y por otro al reflejar el infrarrojo, lo que posibilita, por lo tanto, muchos factores que influyen en la variación estacional y a lo largo de la serie de tiempo. Sin embargo, cuando se toman en cuenta escalas tanto regionales, como ecorregiones, los principales factores que influyen son el tamaño de píxel, la humedad de la vegetación y los diferentes tipos de suelo. En el caso de coberturas antrópicas, tiene que ver el manejo, donde se incluye el riego. En la subregión bsp, por ejemplo, que posee un área de 613 808 ha, se encuentra una combinación de todos los factores mencionados. Para citar algunos, se encuentra el distrito de riego denominado “Arenal-Tempisque” que distribuye aguas para riego y hoy en día cubre 28 000 ha (SENARA, 2018) de la región bsp. Otro factor que incide en la respuesta de las variables biofísicas de esta ecorregión es el hecho de que aquí se encuentran los principales centros urbanos de la provincia de Guanacaste, a saber: Cañas, Santa Cruz, Nicoya y Liberia. El tomar en cuenta que las ecorregiones tienen diferente dinámica antropogénica, es sumamente relevante a la hora de compararlas. Por ejemplo, Fang et al., (2005b) usó fPAR junto con EVI (Enhanced Vegetation Index, por sus siglas en inglés), NDVI y LAI para evaluar la recuperación de la vegetación ante perturbaciones creadas por el ser humano. En la ecorregión bsp, por ejemplo, se da la gran mayoría de los incendios forestales (CONIFOR, 2014), elemento que afecta los procesos biofísicos de la vegetación y por lo tanto los valores de fPAR, NDVI y LAI. Muchos de los incendios forestales y quemas agrícolas están asociados a la preparación del terreno para agricultura y renovación de pastos (CONIFOR, 2014).

En atención a las similitudes de las series de tiempo, todas las ecorregiones poseen tendencias diferentes en magnitud y en forma, la principal diferencia es la magnitud. Tanto bhe2, bmt y bsp presentan un pico mínimo en 2010 y 2011. Por otra parte, bhe1, bhe2 y bsp poseen tendencia a un mínimo en 2010-2011 y máximos en el periodo 2015-2016. Esto es un patrón que se repite en todas las ecorregiones, es decir, que no necesariamente un descenso brusco en la precipitación repercute en una caída de fPAR, LAI o NDVI. Resalta la diferencia entre bhe1 y bhe2, tanto en magnitud como en forma; sin embargo, en algunos sectores de la curva, se comportan en forma muy similar. Esto es importante ya que son definidas como una misma ecorregión (bosques húmedos estacionales) por lo que se supone deberían tener un comportamiento similar. La información presentada en este trabajo puede servir para definir mejor los límites de ecorregiones o bien aportar elementos para una clasificación de los ecosistemas, donde se tome en cuenta el comportamiento de la vegetación de una región.

Para el periodo 2014-2016, la precipitación, para todas las ecorregiones, presenta una tendencia a producir un mínimo (que coincide con un periodo NIÑO), por el contrario, en fPAR y LAI para las ecorregiones bhe1, bhe2 y bsp tienden a un máximo. Esta misma situación se presenta en bhe2 y bsp para el periodo 2004 y 2005. Esto podría indicar que, al menos para la tendencia de fPAR, la propensión de la precipitación no es un factor que se correlaciona directamente con fPAR y LAI a lo largo de toda la curva, es decir, hay periodos de desfase. Sin embargo, esta última aseveración debe corroborarse con un análisis específico, el cual se sale de los alcances de este estudio. Por otra parte, Stöckli et al., (2011), encontraron que la fenología estudiada por medio fPAR y LAI, en ambientes tropicales, no es sensible a condiciones de humedad y variaciones de luz. Sin embargo, el elemento humedad puede ser un aspecto a tomar en cuenta en la ecorregión bmt, donde la precipitación horizontal puede aportar grandes cantidades de agua que bien pueden no ser cuantificadas por productos que brindan información sobre precipitación, como el presentado en este trabajo. Lo anterior podría explicar por qué aunque hay una clara tendencia a disminuir en precipitación NDVI, se inclinan a un máximo en esta ecorregión para el periodo 2015-2016.

Si se comparan fPAR y LAI por ecorregión, estas se comportan muy similarmente en cuanto a su tendencia y estacionalidad, esto indica que no es necesario monitorear ambas variables. Caso contrario lo presenta el NDVI que sí difiere mucho de fPAR y LAI.

La relación entre NDVI y LAI y fPAR, ha sido bien establecida a través de estudios teóricos y empíricos, como por ejemplo en Wang et al., (2005). Sin embargo, la comparación en forma de series de tiempo y por unidades específicas, como en este estudio por ecorregiones, ha sido menos estudiada. La estrategia de usarlas como unidades de comparación biofísicas permitirá realizar el monitoreo no solo de estas variables en el tiempo si no de su relación con otras, como por ejemplo el fenómeno del ENOS, al ser este último punto un objetivo futuro de investigación.

Estudios locales de series de tiempo de variables biofísicas, particularmente aquellos enfocados al análisis de la tendencia y estacionalidad, son escasos en la bibliografía. El estandarizar un método de análisis abre muchas posibilidades para el monitoreo por ecorregiones, ya que la respuesta fenológica es cada vez más relevante para abordar los problemas de productividad ecosistémica. Sin embargo, aún existen desafíos con respecto a las escalas de observación y cuáles son las regiones geográficas de análisis que permitan el modelado del comportamiento fenológico, para permitir pronósticos ecológicos a la luz del cambio climático futuro.

En el método de descomposición de las series de tiempo, las iteraciones de “LOESS” están reguladas por parámetros predefinidos (t.window, s.window y robust), que permiten al modelo tener en cuenta varios factores como la robustez y suavizado, así como controlar el efecto estacional para resaltar la tendencia, lo que ofrece mucha flexibilidad. Lo anterior produce resultados robustos que proporcionan un buen ajuste a los datos, mantiene la transparencia a través de la descomposición de series de tiempo ya que permite examinar las diferentes curvas, al resaltar en la tendencia de cada una de las series. Por otra parte, la fácil integración de la novedosa plataforma de Google Earth Engine y el lenguaje de programación R, con paquetes y funciones, facilitan la obtención y posterior descomposición de las series de tiempo, con lo cual se simplifican las labores de monitoreo.

5.Conclusiones

Este estudio presenta dos formas de análisis de variables biofísicas, la primera es la comparación local de las mismas por ecorregiones y la segunda, el análisis de estas variables en forma de series de tiempo, que de forma descompuesta, facilitan el análisis de la tendencia. Pese a que se han desarrollado muchos estudios sobre variables biofísicas, estos son regionales y se necesitan más estudios locales para mejorar la comprensión del efecto en cuestiones como el uso y cobertura de la tierra en el pasado, el manejo de los cultivos, el efecto del riego y las formas de producción actuales. Las variables climáticas que se relacionan con la producción primaria por ecorregión, facilitan el monitoreo y posterior modelado biogeográfico y ecológico en cuatro ecorregiones de la provincia de Guanacaste. Es necesario, además, buscar la correlación de la tendencia de las variables biofísicas estudiadas, con, por ejemplo, el efecto ENOS.

Las técnicas de sensoramiento remoto presentan muchas ventajas en la estimación de LAI, fPAR y NDVI sobre los métodos tradicionales de medición de campo y proporcionan el potencial para analizar las variaciones temporales y espaciales de estas variables a diferentes escalas y en unidades espaciales específicas como ecorregiones.

La región estudiada se caracteriza por su adaptación a la fuerte estacionalidad. Sin embargo, los estudios de modelización ecológica en Guanacaste no pueden obviar la dinámica del uso y cobertura de la tierra de las diferentes ecorregiones. Por ejemplo, la bsp es la región más dinámica en cuanto a producción agropecuaria y ganadera, esto hace que el monitoreo del bioclima deba tomar en cuenta estos aspectos.

El uso de los productos globales de LAI, fPAR y NDVI existentes está limitado por su resolución espacial (500 m). Ninguno de los productos globales de LAI, fPAR o NDVI ha sido comparado con los datos para diferentes ecosistemas, tanto productivos como naturales.

La investigación futura podría enfocarse, por un lado, en el monitoreo de variables biofísicas por ecorregión, y por el otro, a la correlación de los productos como LAI, fPAR y NDVI con fenómenos como el ENOS u otros índices que influencian el clima de una región como la de Guanacaste. Para esto las novedosas herramientas como Earth Engine y Climate Engine han mostrado ser elementos invaluables que facilitan las labores de obtención de las series de tiempo.

Las ecorregiones analizadas en Guanacaste presentan matrices heterogéneas de usos del suelo, de ecosistemas naturales, tipos de suelo y topografía, por lo que evaluar el historial de su tendencia y estacionalidad a partir de variables biofísicas tiene sentido para evaluar los cambios en productividad ecosistémica.

Dada la similaridad de las curvas de fPAR y LAI en tendencia y estacionalidad, los autores sugieren no evaluar ambas variables para estudios futuros.

Finalmente, análisis locales de series de tiempo de variables biofísicas con metodologías de análisis estandarizados abren muchas posibilidades para el monitoreo por ecorregiones, ya que la respuesta fenológica es cada vez más relevante para abordar los problemas de productividad ecosistémica. Sin embargo, aún existen desafíos con respecto a las escalas de observación y la definición de cuáles son las regiones de análisis que permitan el modelado del comportamiento fenológicas para consentir pronósticos ecológicos a la luz del cambio climático futuro.

6.Agradecimientos

Los autores agradecen profundamente el apoyo brindado por el Consejo Nacional para Investigaciones Científicas y Tecnológicas (CONICIT) y al Ministerio de Ciencia Tecnología y Telecomunicaciones (MICITT) al financiar el proyecto “Seguridad y Sustentabilidad Hídrica en Guanacaste: Innovaciones en el Análisis Hidroclimático Tradicional y su Aplicación en el Manejo Integral del Recurso Hídrico en un Clima Cambiante”. Finalmente, se agradece a la Revista y a las personas revisoras anónimas por las observaciones realizadas al manuscrito.

7.Referencias

Amthor, J. y Baldocchi, D. (2001). Terrestrial Higher Plant Respiration and Net Primary Production. Terrestrial Global Productivity, 33-59. https://doi.org/10.1016/b978-012505290-0/50004-1

Asner, G. P., Scurlock, J. M. O., y A. Hicke, J. (2003). Global synthesis of leaf area index observations: implications for ecological and remote sensing studies. Global Ecology and Biogeography, 12(3), 191-205. https://doi.org/10.1046/j.1466-822X.2003.00026.x

Baldocchi, D. D. y Wilson, K. B. (2001). Modeling CO2 y water vapor exchange of a temperate broadleaved forest across hourly to decadal time scales. Ecological Modelling, 142(1), 155-184. https://doi.org/10.1016/s0304-3800(01)00287-3

Bicheron, P. (1999). A Method of Biophysical Parameter Retrieval at Global Scale by Inversion of a Vegetation Reflectance Model. Remote Sensing of Environment, 67, 251-266. https:// doi.org/10.1016/S0034-4257(98)00083-2. https://doi.org/10.1016/s0034-4257(98)00083-2

Blümel, B. y Reimer, E. (2009). Validation of boundary layer parameters of climate model REMO: estimation of leaf area index from NOAA-AVHRR data for the Baltimos region. Theoretical and Applied Climatology, 98, 1-11. https://doi.org/10.1007/s00704-009-0164-3

Bonan, G. (2016). Ecological Climatology. Concepts and Applications. Cambridge University Press, 3rd edition edition. https://doi.org/10.1017/CBO9780511805530

Box, G. E. P. y Pierce, D. A. (1970). Distribution of residual autocorrelations in autoregressive- integrated moving average time series models. Journal of the American Statistical Association, 65(332), 1509-1526. https://doi.org/10.2307/2284333

Brockwell, P. J. y Davis, R. A. (2002). Introduction to time series and forecasting. Springer- Verlag, New York.

Brown, L., Chen, J., Leblanc, S., y Cihlar, J. (2000). Shortwave infrared correction to the simple ratio: An image and model analysis. Remote Sens. Environ, 77, 16-25. https://doi.org/10.1016/S0034-4257(99)00035-8

Brown, M. E., Pinzon, J. E., Didan, K., Morisette, J. T., y Tucker, C. J. (2006). Evaluation of the consistency of long-term NDVI time series derived from AVHRR, Spot-Vegetation, Seawifs, MODIS, and Landsat ETM+ sensors. IEEE Transactions on Geoscience and Remote Sensing, 44(7), 1787-1793. https://doi.org/10.1109/TGRS.2005.860205

Chen, J. (1996). Evaluation of vegetation indices and a modified simple ratio for boreal applications. Canadian Journal of Remote Sensing, 22(3), 229-242. https://doi.org/ 10.1080/07038992.1996.10855178

Chen, J., Liu, J., Cihlar, J., y Goulden, M. (1999). Daily canopy photosynthesis model through temporal and spatial scaling for remote sensing applications. Ecological Modelling, 124(2):99-119. https://doi.org/10.1016/S0304-3800(99)00156-8

Chen, P. Y., Fedosejevs, G., Tiscareño-LóPez, M., y Arnold, J. G. (2006). Assessment of MODIS-EVI, MODIS-NDVI and Vegetation-ndvi composite data using agricultural measurements: An example at corn fields in western Mexico. Environmental Monitoring and Assessment, 119(1), 69-82. DOI: https://doi.org/10.1007/s10661-005-9006-7

Choudhury, B. J. (1987). Relationships between vegetation indices, radiation absorption, and net photosynthesis evaluated by a sensitivity analysis. Remote Sensing of Environment, 22(2), 209-233. https://doi.org/10.1016/0034-4257(87)90059-9

Cleveland, R., Cleveland, W., McRae, J., y Terpenning, I. (1990). A seasonal-trend decomposition procedure based on loess (with discussion). Journal of Official Statistics, 6(1).

CONIFOR. (2014). Estrategia nacional de manejo integral del fuego en Costa Rica 2012-2021. Technical report, Comisión Nacional sobre Incendios Forestales CONIFOR, SINAC.

Davis, J. (2003). Statistics and Data analysis in Geology, John Wiley and Sons.

Didan, K. (2015). Mod13q1 modis/terra vegetation indices 16-day l3 global 250m sin grid v006.. Technical report, NASA EOSDIS Land Processes DAAC. https://doi.org/10.5067/modis/ mod13q1.006.

Eklundh, L., Harrie, L., y Kuusk, A. (2001). Investigating relationships between lysat etm+ sensor data and leaf area index in a boreal conifer forest. Remote Sensing of Environment. https://doi.org/10.1016/S0034-4257(01)00222-X

Fang, H., Liang, S., McClaran, M., Van Leeuwen, W., Drake, S., Marsh, S., Thomson, A., Izaurralde, R., y Rosenberg, N. (2005a). Biophysical characterization and management effects on semiarid rangeland observed from landsat etm+ data. IEEE Transactions on Geoscience and Remote Sensing, 43(1):125-133. https://doi.org/10.1109/tgrs.2004.839813

Fang, H., Liang, S., McClaran, M. P., van Leeuwen, W. J. D., Drake, S., Marsh, S. E., Thomson, A. M., Izaurralde, R. C., y Rosenberg, N. J. (2005b). Biophysical characterization and management effects on semiarid rangeland observed from Landsat ETM+ data. IEEE Transactions on Geoscience and Remote Sensing, 43(1), 125-134. https://doi.org/10.1109/tgrs.2004.839813

Feng, D., Chen, J., Plummer, S., Chen, M., y Pisek, J. (2006). Algorithm for global leaf area index retrieval using satellite imagery. IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, 44(8). https://doi.org/10.1109/TGRS.2006.872100

Fensholt, R., Sandholt, I., y Schultz Rasmussen, M. (2004). Evaluation of modis lai, fapar and the relation between fapar and ndvi in a semi-arid environment using in situ measurements. Remote Sensing of Environment, 91, 490-507. https://doi.org/10.1016/j.rse.2004.04.009

Forkel, M., Carvalhais, N., Verbesselt, J., Mahecha, M. D., Neigh, C. S., y Reichstein, M. (2013). Trend change detection in ndvi time series: Effects of inter-annual variability and methodology. Remote Sensing, 5(5), 2113-2144. https://doi.org/10.3390/rs5052113

Friedl, M. A., Davis, F. W. Michaelsen, J., y Moritz, M. A. (1995). Scaling and uncertainty in the relationship between the ndvi and land surface biophysical variables: an analysis using a scene simulation model and data from fife. Remote Sensing of Environment, 54, 233- 246. https://doi.org/10.1016/0034-4257(95)00156-5

Fuller, W. A. (1996). Introduction to Statistical Time Series. John Wiley & Sons, Inc. http://dx.doi.org/10.1002/9780470316917

Funk, C., Peterson, P., Landsfeld, M.F., Pedreros, D., Verdin, J.P., Rowland, J., Romero, B., Husak, G., Michaelsen, J., y Verdin, A. (2014). A quasi-global precipitation time series for drought monitoring: U.S. Geological Survey Data Series 832, p. 4 http://dx.doi.org/10.3133/ds832

Garrigues, S., Lacaze, R., Baret, F., Morisette, J. T., Weiss, M., Nickeson, J. E., Fernandes, R., Plummer, S., Shabanov, N. V., Myneni, R. B., Knyazikhin, Y., y Yang, W. (2008). Validation and intercomparison of global leaf area index products derived from remote sensing data. Journal of Geophysical Research: Biogeosciences, 113(G2), n/a–n/a. G02028. https://doi.org/10.1029/2007JG000635

Gobron, N. y Verstraete, M. M. (2009). Fraction of absorbed photosynthetically active radiation (fapar): Assessment of the status of the development of the standards for the terrestrial essential climate variables. Technical Report Version 8, Global Terrestrial Observing System-FAO.

Gorelick, N., Hancher, M., Dixon, M., Ilyushchenko, S., Thau, D., y Moore, R. (2017). Google Earth Engine: Planetary-scale Geospatial Analysis for Everyone. Remote Sensing of Environment, 202 (1), 18-27. https://doi.org/10.1016/j.rse.2017.06.031

Gower, S. T., Krankina, O., Olson, R. J., Apps, M., Linder, S., y Wang, C. (2001). Net primary production and carbon allocation ppattern of boreal forest ecosystems. Ecological Applications, 11(5), 1395–1411. https://doi.org/10.1016/j.rse.2017.06.031

Hobbs, T. J. (1995). The use of NOAA-AVHRR NDVI data to assess herbage production in the arid rangelands of Central Australia. International Journal of Remote Sensing, 16(7), 1289- 1302. https://doi.org/10.1080/01431169508954477

Hobijn, B., Franses, P. H., y Ooms, M. (1998). Generalizations of the kpss-test for stationarity. Statistica Neerlandica, 58 (4), 483-502. https://doi.org/10.1111/j.1467-9574.2004.00272.x

Huete, A., Didan, K., Miura, T., Rodriguez, E., Gao, X., y Ferreira, L. (2002). Overview of the Radiometric and Biophysical Performance of the MODIS Vegetation Indices. Remote Sensing of Environment, 83(1-2):195–213. https://doi.org/10.1016/S0034-4257(02)00096-2

Huete, A., Didan, R., van Leeuwen, W., Miura, T., y Glenn, E. (2008). Modis vegetation indices. in: Land remote sensing and global environmental change: Nasa’s earth observing system and the science of aster and modis. Land Remote Sensing and Global Environmental Change: NASA’s Earth Observing System and the Science of ASTER and MODIS. https://doi.org/10.1007/978-1-4419-6749-7_26

Hunt, E. R., Kelly, R. D., Smith, W. K., Fahnestock, J. T., Welker, J. M., y Reiners, W. A. (2004). Estimation of carbon sequestration by combining remote sensing and net ecosystemexchange data for northern mixed-grass prairie and sagebrush–steppe ecosystems. Environmental Management, 33(1): S432-S441. https://doi.org/10.1007/s00267-003-9151-0

Huntington, J., Hegewisch, K., Daudert, B., Morton, C., Abatzoglou, J., McEvoy, D., y Erickson, T. (2017). Climate Engine: Cloud Computing of Climate and Remote Sensing Data for Advanced Natural Resource Monitoring and Process Understanding. Disponible en: https://clim-engine.appspot.com/.X. https://doi.org/10.1175/bams-d-15-00324.1

Hyndman, R. J. y Athanasopoulos, G. (2013). Forecasting: principles and practice. Texts. Disponible en: https://otexts.com/fpp2/

Hyndman, R. J. y Khandakar, Y. (2008). Automatic time series forecasting: the forecast package for R. Journal of Statistical Software, 26(3), 1-22.

Iwata, H., Ueyama, M., Iwama, C., y Harazono, Y. (2013). Variations in fraction of absorbed photosynthetically active radiation and comparisons with modis data in burned black spruce forests of interior Alaska. Polar Science, 7(2):113 - 124. https://doi.org/10.1016/j.polar.2013.03.004

Jackson, D. R. y Huete, R. A. (1991). Intepreting vegetation indices. Preventive Veterinary Medice, 11:185-200. https://doi.org/10.1016/s0167-5877(05)80004-2

Jacquemoud, S., Verhoef, W., Baret, F., Bacour, C., Zarco-Tejada, P., Asner, G., François, C., Ustin, S. (2009). Use of prospect+sail to estimate canopy biochemistry at different scales. Remote Sensing of Environment, 113(1), S56-S66. https://doi.org/10.1016/j.rse.2008.01.026

Jacquemoud, S. y Ustin, S. (2003). Application of radiative transfer modes to moisture content estimation and burned land mapping. In Proceedings of the 4th International Workshop on Remote Sensing and GIS applications to Forest Fire management, Innovative concepts and methods in fire danger estimation.

Jonckheere, I., Fleck, S., Nackaerts, K., Muys, B., Coppin, P., Weiss, M., y Baret, F. (2004). Review of methods for in situ leaf area index determination part i. theories, sensors and hemispherical photography. Agricultural and Forest Meteorology, 121(1-2), 19-35. https://doi.org/10.1016/j.agrformet.2003.08.027

Kappas, M. W. y Propastin, P. A. (2012). Review of available products of leaf area index and their suitability over the formerly Soviet Central Asia. Journal of Sensors, 2012. https://doi.org/10.1155/2012/582159

Knyazikhin, Y., Martonchik, J. V., Diner, D. J., Myneni, R. B., Verstraete, M., Pinty, B., y Gobron, N. (1998). Estimation of vegetation canopy leaf area index and fraction of absorbed photosynthetically active radiation from atmosphere-corrected misr data. Journal of Geophysical Research: Atmospheres, 103(D24), 32239-32256. https://doi.org/10.1029/98jd02461

Künzer, C., Dech, S., y Wagner, W. (2015). Remote Sensing Time Series: Revealing Land Surface Dynamics. 22. Springer International Publishing. https://doi.org/10.1007/978-3-319-15967-6

Lotsch, A., Tian, Y., Friedl, M. A., y Myneni, R. B. (2003). Land cover mapping in support of LAI and fpar retrievals from eos-modis and misr: Classification methods and sensitivities to errors. International Journal of Remote Sensing, 24(10), 1997-2016. https://doi.org/10.1080/01431160210154858.X

Luyssaert, S., Inglima, I., Jung, M., Richardson, A., Reichstein, M., Papale, D., Piao, S., Schulze, E., Wingate, L., Matteucci, G., Aragao, L., Aubinet, M., Beers, C., Bernhofer, C., Black, G., Bonal, D., Bonnefond, J., Chambers, J., Ciais, P., Cook, B., Davis, K., Dolman, A., Gielen, B., Goulden, M., Grace, J., Granier, A., Grelle, A., Griffis, T., Grünwald, T., Guidolotti, G., Hanson, P., Harding, R., Hollinger, D., Hutyra, L.R., Kolari, P., Kruijt, B., Kutsch, W., Lagergren, F., Laurila, T., Law, B., Le Maire, G., Lindroth, A., Loustau, D., Malhi, Y., Mateu, J., Migliavacca, M., Misson, L., Montagnani, L., Moncrieff, J., Moors, E., Munger, J., Nikinmaa, E., Ollinger, S., Pita, G., Rebmann, C., y Roupsard, O. (2007). CO2 Balance of Boreal, Temperate, and Tropical Forests Derived from a Global Database. Global Change Biology, 13(12), 2509-2537. https://doi.org/10.1111/j.1365-2486.2007.01439.x

min WANG, F., feng HUANG, J., lin TANG, Y., y zhen WANG, X. (2007). New vegetation index and its application in estimating leaf area index of rice. Rice Science, 14(3), 195- 203. https://doi.org/10.1016/s1672-6308(07)60027-4

Morisette, J. T., Richardson, A. D., Knapp, A. K., Fisher, J. I., Graham, E. A., Abatzoglou, J., Wilson, B. E., Breshears, D. D., Henebry, G. M., Hanes, J. M., y Liang, L. (2009). Tracking the Rhythm of the Seasons in the Face of Global Change: Phenological Research in the 21st Century. Frontiers in Ecology and the Environment, 7(5), 253-260. https://doi.org/10.1890/070217

Myneni, R. (2015). Mod15a2h modis/terra leaf area index/fpar 8-day l4 global 500m sin grid v006. Technical report, NASA EOSDIS Land Processes DAAC. https://doi.org/10.5067/modis/mod15a2h.006.

Myneni, R., Hoffman, S., Knyazikhin, Y., Privette, J., Glassy, J., Tian, Y., Wang, Y., Song, X., Zhang, Y., Smith, G., Lotsch, A., Friedl, M., Morisette, J., Votava, P., Nemani, R., y Running, S. (2001). Global products of vegetation leaf area and fraction absorbed par from year one of modis data. Remote Sensing of Environment, 83, 214-231.

Myneni, R., Maggion, S., Iaquinta, J., Privette, J., Gobron, N., Pinty, B., Kimes, D., Verstraete, M., y Williams, D. (1995). Optical remote sensing of vegetation: Modeling, caveats, and algorithms. Remote Sensing of Environment, 51(1), 169-188. Remote Sensing of Land Surface for Studies of Global Chage. https://doi.org/10.1016/0034-4257(94)00073-v

Myneni, R., Ramakrishna, R., Ramakrishna, N., y Running, S. (1997). Estimation of global leaf area index and absorbed par using radiative transfer models. IEEE Transactions on Geoscience and Remote Sensing 35(6):1380-1393. https://doi.org/10.1109/36.649788

Olson, D., & Dinerstein, E. (2002). The Global 200: Priority Ecoregions for Global Conservation. Annals of the Missouri Botanical Garden, 89(2), 199-224. https://doi.org/10.2307/3298564

Phillips, P. C. B. y Perron, P. (1988). Testing for a unit root in time series regression. Biometrika, 75(2), 335-346. https://doi.org/10.1093/biomet/75.2.335

Price, J. (1994). How unique are spectral signatures? Remote Sensing of Environment, 49(3), 181- 186. https://doi.org/10.1016/0034-4257(94)90013-2

Programa REDD/CCAD-GIZ-SINAC. 2015. Inventario Nacional Forestal de Costa Rica 2014-2015. Resultados y Caracterización de los Recursos Forestales. Preparado por: Emanuelli, P., Milla, F., Duarte, E., Emanuelli, J., Jiménez, A. y Chavarría, M.I. Programa Reducción de Emisiones por Deforestación y Degradación Forestal en Centroamérica y la República Dominicana (REDD/CCAD/GIZ) y Sistema Nacional de Áreas de Conservación (SINAC) Costa Rica. San José, Costa Rica, 380 p.

Propastin, P. A. (2009). Spatial non-stationarity and scale-dependency of prediction accuracy in the remote estimation of lai over a tropical rainforest in Sulawesi, Indonesia. Remote Sensing of Environment, 113(10), 2234-2242. https://doi.org/10.1016/j.rse.2009.06.007

R Core Team. (2018). R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria.

Running, S. W., Nemani, R. R., Heinsch, F. A., Zhao, M., Reeves, M., y Hashimoto, H. (2004). A continuous satellite-derived measure of global terrestrial primary production. BioScience, 54(6), 547-560. https://doi.org/10.1641/0006-3568(2004)054[0547:acsmog]2.0.co;2

SENARA. (2018). ¿Qué es el DRAT? Disponible en: http://www.senara.or.cr/drat/.X

Shumway, R. y Stoffer, D. (2017). Time Series Analysis: Using the R Statistical Package. Free Dog Publishing.

Stöckli, R. T., Rutishauser, I., Liniger, B. M., y Denning, A. (2011). A global reanalysis of vegetation phenology, J. Geophys. Res., 116(G03020). https://doi.org/10.1029/2010JG001545.

Trenberth, K. E. (1997). The definition of El Niño. Bulletin of the American Meteorological Society, 78(12), 2771-2777. https://doi.org/10.1175/1520-0477(1997)078<2771:tdoeno>2.0.co;2

Tucker, C. J. (1979). Red and photographic infrared linear combinations for monitoring vegetation. Remote Sensing of Environment, 8(2), 127-150. https://doi.org/10.1016/0034-4257(79)90013-0

Varpe, Ø. (2017). Life history adaptations to seasonality. Integrative and Comparative Biology, 57(5), 943-960. https://doi.org/10.1093/icb/icx123.

Viña, A. y Gitelson, A. (2005). New developments in the remote estimation of the fraction of absorbed photosynthetically active radiation in crops. Geophysical Research Letters, 32(17), 1-4. https://doi.org/10.1029/2005GL023647

Wang, Q., Adiku, S., Tenhunen, J., y Granier, A. (2005). On the relationship of ndvi with leaf area index in a deciduous forest site. Remote Sensing of Environment, 94(2), 244-255. https://doi.org/10.1016/j.rse.2004.10.006

8.Apéndices

Apéndice 1. Amplitud de NDVI para las cuatro ecorregiones en Guanacaste. Se aprecia muy bien la diferencia entre los bosques secos del pacífico (bsp) y los bosques húmedos estacionales (bhe1 y bhe2). Sin embargo, datos de amplitud en NDVI podrían servir para afinar los límites de las ecorregiones.

Apéndice 2. Imagen de fase de NDVI para las cuatro ecorregiones en Guanacaste. Se aprecia muy bien la diferencia entre los bosques montanos de Talamanca (bmt) y el resto de ecorregiones. El usar la fase del NDVI puede ser un elemento a tomar en cuenta para futuras procesos de clasificación de ecosistemas.

Apéndice 3. Distribución de la precipitación total mensual en mm durante el periodo 1981-2017 para las cuatro ecorregiones en Guanacaste. Datos del producto “CHIRPS” basados en (Funk et al., 2014) y obtenidos de la aplicación Climate Engine (Huntington et al., 2017).

Con esta misma serie de tiempo se obtuvieron los datos mensuales para todos los años, serie que se presenta descompuesta en estacionalidad, residuos y tendencia en el Apéndice 4. Para las ecorregiones bhe2, bmt y bsp entre el 2010 y 2011 se observa un pico elevado de precipitación que coincide con un evento NIÑA para esos años. En las cuatro ecorregiones hay un pico mínimo entre 2015 y 2016 que coincide con un evento NIÑO muy fuerte.

Apéndice 4. Serie de tiempo (periodo 2002-2017) con estacionalidad, tendencia y componentes irregulares para la precipitación en las cuatro ecorregiones en Guanacaste. Datos del producto “CHIRPS” basados en (Funk et al., 2014) y obtenidos de la aplicación Climate Engine (Huntington et al., 2017).

Apéndice 5. Series de tiempo mensual (periodo 2002-2017) para la precipitación en mm en las cuatro ecorregiones en Guanacaste. Datos del producto “CHIRPS” basados en (Funk et al., 2014) y obtenidos de la aplicación Climate Engine (Huntington et al., 2017).

Apéndice 6. Líneas de tendencia de precipitación, fPAR, NDVI y LAI para las cuatro ecorregiones, en 2002-2017 en las cuatro ecorregiones en Guanacaste.

El Niño y la Oscilación del sur, también conocido como ENOS es una fluctuación periódica en la temperatura de la superficie del mar (El Niño) y la presión atmosférica de la atmósfera que lo cubre (Oscilación del sur) a través del océano Pacífico ecuatorial (Trenberth, 1997).

Apéndice 7. Serie de tiempo de la Anomalía de la temperatura de la superficie del mar (SST, por sus siglas en inglés) en la región NIÑO 3 (1990-2017).

---