El primer paso de la metodología aplicada en este trabajo consistió en la selección de una serie de incendios representativos de las comunidades vegetales afectadas por el fuego en Aragón que, además, reunieran ciertos requisitos en cuanto a sus dimensiones (superficie quemada > 500 ha) y distribución espacial (ubicación dentro de un único tile de MODIS). Para ello se contó con la base de datos de incendios forestales del Gobierno de Aragón proporcionada por el Servicio de Gestión de los Incendios Forestales y Coordinación (Dirección General de Gestión Forestal) y de las cartografías en formato digital generadas en el contexto del proyecto GA-LC-042/2011 (Caixa-DGA), mediante técnicas de teledetección y fuentes de datos oficiales.

La selección de los incendios se centró en aquellos ocurridos entre 2000 y 2010, con el propósito de realizar un seguimiento anual de los valores de LAI/fPAR y del impacto del fuego en los primeros meses. Además, la selección se completó incluyéndose incendios antiguos, desde 1975, cuya localización coincidiera con el tile 17h y 4v, coordenadas del extracto de los productos MODIS, con el objeto de evaluar el estado actual de los valores LAI/fPAR en incendios ocurridos en diferentes periodos de tiempo.

En total, se seleccionaron 18 perímetros, espacialmente repartidos en las provincias de Zaragoza y Huesca (Aragón) (Tabla 1 y Figura 1). Debido a la resolución del producto MCD15A2 de MODIS (1 km²), se seleccionaron sólo los píxeles incluidos completamente en el perímetro del incendio, eliminándose los de borde o todos aquellos que, aun estando dentro del perímetro de la superficie afectada, presentaran sectores con zonas de vegetación sin afectar.

La Figura 1 recoge la distribución espacial de los incendios forestales analizados. Se trata de incendios que se distribuyen a lo largo de un eje NE-SO en la zona occidental de Aragón, pudiéndose agrupar en 4 grandes marcos geográficos: incendios localizados en la Ibérica zaragozana, incendios del CENAD (Centro Nacional de Adiestramiento “San Gregorio”), incendios de las sierras prepirenaicas, e incendios de Zuera y los Montes de Castejón. Estos marcos de referencia espacial presentan cierta homogeneidad en relación a las formaciones vegetales afectadas, condiciones fisiográficas y dominios climáticos.

Los incendios con los códigos 6, 7, 11, 13, 17, 19, 26, 29 y 31 recogen amplios sectores quemados integrados en las sierras occidentales del Prepirineo oscense y zaragozano, cubiertos por comunidades de carácter submediterráno dominadas por Quercus gr. cerrioides (Willk y Costa) y, en menor medida, por Quercus ilex L.; pinares de Pinus sylvestris L. y Pinus nigra (Arnold) y formaciones arbustivas dominadas por Buxus sempervirens L., Genista scorpius L., etc. (Pérez-Cabello, 2002; Ruiz de la Torre, 1990a). Los incendios con los códigos 10 y 34 recogen información de los incendios de los Montes de Zuera y Castejón de Valdejasa, siendo representativos de formaciones arbóreas dominadas por Pinus halepensis L., acompañadas por Q. coccifera L. y otras especies típicas de matorral mediterráneo. Los incendios con códigos 32 y 38 se localizan en el CENAD y remiten a formaciones vegetales en su mayoría de carácter subarbustivo y pastizales xerofíticos.

Finalmente, los incendios con los códigos 23, 24, 25, 36 y 37 se localizan en diferentes sectores de las estribaciones meridionales del macizo del Moncayo y otras sierras zaragozanas. En este caso, las zonas afectadas estaban cubiertas de manera predominante por matorrales incluidos en la alianza mediterránea Rosmarino-Ericion, encinares y plantaciones de Pinus pinaster (Aiton) y Pinus halepensis L. (Ruiz de la Torre, 1990b; Villar, 1990).

Desde el punto de vista de la estrategia reproductiva tras el fuego, esta amplia muestra de ecosistemas contiene especies vegetales con tres tipos básicos de respuesta al fuego (Pausas, 2012): (a) Reclutadores obligados dominados por P. halepensis, especie dotada de un extenso banco de semillas aéreo (serotinia) que posibilita una recolonización de tipo endógeno; (b) Sistemas pasivos de defensa, en este caso se trata de comunidades dominadas por P. nigra y P. sylvestris que mantienen una estrategia resistente al fuego mediante sistemas pasivos y germinación exógena; y (c) Rebrotadores obligados, reúnen a fanerófitos y macro-fanerófitos pertenecientes al género Quercus (Q. ilex, Q. gr. cerrioides, Q. coccifera) con estrategia de rebrote mediante yemas adventicias.

Una vez recabados todos los perímetros de los incendios de interés se procedió a la descarga de los datos MODIS a través del servidor Earth Explorer del USGS (http://earthexplorer.usgs.gov/). Tal y como se cita más arriba, al descargar los datos se seleccionó el tile o imagen de 17h y 4v de la rejilla de la proyección sinusoidal de MODIS.

Se descargaron dos compuestos estacionales de 8 días, primera semana de mayo y septiembre, por considerarse claves en el registro del contraste entre la actividad fotosintética de la vegetación a finales de primavera y el descenso de la misma al final del verano, justo antes de que se generalicen las lluvias otoñales (López et al., 2007). Además, en septiembre la geometría de iluminación presenta menos problemas que durante los meses siguientes como consecuencia del descenso del ángulo de elevación solar.

Por otro lado, se utilizaron los parámetros de calidad (Quality flags –QC–) del propio producto con el propósito de eliminar posibles píxeles alterados por la presencia de nubes, sombras de nube, etc. A partir de aquí se tomarán las siguientes denominaciones para los compuestos y variables: LAI_sep y LAI_may para los compuestos de septiembre y mayo, respectivamente, del índice de área foliar; y fPAR_sep y fPAR_may para los compuestos de la variable fracción de la radiación fotosintéticamente activa. Los valores que adquieren estos dos parámetros pueden oscilar entre 0 y 8 para el caso del LAI y entre 0 y 1 para fPAR.

La diferente información recopilada a partir de los productos MODIS se analiza en dos vertientes diferentes: (1) análisis estático de los valores de LAI/fPAR correspondientes a zonas quemadas a lo largo del periodo 1975-2010; (2) seguimiento multianual en 6 grandes incendios ocurridos en la década de 2000, representativos de condiciones ambientales diferentes.

(1) El análisis estático de los valores de LAI/fPAR (compuestos de 2010) referidos a zonas afectadas en diferentes años a lo largo del periodo (1976-2010), nos permite analizar el factor tiempo en la recuperación de los valores de LAI/fPAR. En los compuestos de 2010 quedan recogidos todos los incendios del periodo elegido.

Debido al hecho de que no se dispone de zonas quemadas en todos y cada uno de los años en este periodo, los diferentes incendios han sido agrupados en intervalos temporales en función del año que ocurrieron. A continuación se enumeran dichas agrupaciones: I₁ = superficies quemadas con un año de antigüedad en relación con el último compuesto disponible (2010); I₂ = superficies quemadas con dos años de antigüedad; I₁₀ superficies quemadas entre 2 y 10 años de antigüedad; I₂₀ = superficies quemadas entre 10 y 20 años de antigüedad; I₃₀ = superficies quemadas entre 20 y 30 años de antigüedad; I₃₅ = superficies quemadas con más de 35 años de antigüedad.

El perfil asimétrico de los intervalos se relaciona con el ajuste exponencial entre tiempo y NDVI descrito en muchos tipos de vegetación de tipo mediterránea (Viedma et al., 1997), caracterizado por un rápido incremento en los primeros años tras la perturbación seguido de etapas con un crecimiento más moderado, hasta llegar al comportamiento asintótico a medio plazo.

(2) El segundo de los análisis se trata de un estudio a mayor resolución temporal (análisis multianual) considerando de manera individual 6 grandes incendios producidos en la década de 2000 (a partir de 2002 se dispone de cartografía del producto MODIS MCD13A2). Este segundo análisis permite comprobar los efectos inmediatos del fuego sobre los valores de LAI/fPAR y caracterizar la respuesta de la vegetación, a una resolución temporal de mayor detalle, en zonas quemadas inscritas en diferentes marcos de referencia espacial de la región aragonesa.

Desde el punto de vista de los análisis estadísticos aplicados, mediante análisis de varianza (ANOVA) y la prueba post hoc de Tamhane T2, con gráficos de intervalos de confianza para la media, se comprueba si las medias de los valores de LAI/fPAR presentan diferencias estadísticamente significativas entre agrupaciones temporales. Además se realizan análisis de regresión simple, utilizando el coeficiente de determinación r² como medida de ajuste, para identificar tendencias en el tiempo y el grado de ajuste entre los valores de fPAR y LAI.

En relación con los niveles de fPAR/LAI recogidos en zonas quemadas de Aragón en el periodo comprendido entre 1975 y 2010, en la Tabla 2 se recogen estadísticos descriptivos referentes a los compuestos de mayo y septiembre (fPAR_may/sep y LAI_may/sep) de 2010. Dichos valores han sido extraídos en diferentes píxeles en función de la fecha en la que se produjo el fuego (desde hace un año, hasta píxeles que se vieron afectados por el fuego hace más de 35 años). En total la información analizada proviene de un total 127 píxeles que se corresponden con 127 km² de área quemada procedentes de los 18 incendios considerados en este trabajo.

En lo relativo a los intervalos temporales que se han propuesto para el análisis conjunto de las superficies quemadas, los valores más bajos, tanto de LAI como de fPAR, los encontramos lógicamente en el intervalo I₁, un año después del incendio, en los compuestos de septiembre (fPAR_sep = 0,16; LAI_sep = 0,24). Estos valores, por tanto, pueden ser representativos del grado de regeneración de las estructuras vegetales en el primer año tras el fuego; evidentemente, se trata de un valor medio proveniente de varias zonas quemadas sin tener en cuenta diferentes factores de carácter coyuntural y/o estructural, tales como la severidad o el tipo de formación vegetal afectada, respectivamente.

En el otro extremo, los valores más elevados (fPAR = 0,50, en cualquiera de los compuestos; LAI = 1,13/1,40 para septiembre y mayo respectivamente) se corresponden con las superficies quemadas con mayor antigüedad (I₃₅). Aunque los valores en los intervalos extremos no son estrictamente comparables, dado que agrupan incendios diferentes en cuanto a las formaciones vegetales afectadas, contexto ambiental, etc., resulta significativo el hecho de que, 35 años tras el fuego, podamos encontrar valores que en el caso de fPAR triplican y en el caso de LAI quintuplican a los recogidos un año tras el fuego.

La mayor cuantía de los valores de LAI/fPAR en I_30-35 puede interpretarse en términos de recuperación de la productividad, incremento del nivel de intercambio de los flujos de masas y energía en el sistema atmósfera-biosfera y, desde el punto de vista de la protección del suelo, mejora en los parámetros controladores de la interceptación de la lluvia. No obstante, estos valores deben ser interpretados con cierta cautela ya que numerosos trabajos de verificación realizados hasta la fecha sobre los productos MODIS alertan de la sobreestimación de los registros del producto MCD15A2 en comparación con observaciones de campo (Steinberg et al., 2006). Durá et al. (2012) identificaron una sobreestimación del 38% durante el año fenológico 2009-2010 en ambientes de dehesa en el que se producen dinámicas de LAI muy acusadas debido a la contribución de la componente herbácea de estos ecosistemas.

Esta disfunción se relaciona especialmente con los problemas de clasificación de los tipos de cubierta (Serbin et al., 2013) ya que, una mala asignación del bioma dentro de los algoritmos de LAI/fPAR, puede introducir errores superiores al 20% (Myneni et al., 2002). En este sentido, a nivel estacional, se han identificado diferencias estadísticamente significativas (p-value <0,05) en los valores de fPAR y LAI, entre los compuestos de mayo y septiembre, presentándose siempre guarismos más elevados en los compuestos de mayo (1,5 veces mayores en el caso del LAI y 1,2 en el caso de fPAR), al situarse en la fase primaveral del ciclo fenológico de la vegetación, momento en el que el desarrollo y la complejidad de su estructura foliar es mayor.

Por otro lado, hay que reseñar la alta variabilidad que presentan los valores registrados por intervalo temporal debido (1) a la heterogeneidad de las formaciones vegetales afectadas y de sus estrategias reproductivas; (2) a las diferencias en la coyuntura pluvio-térmica experimentada por cada zona quemada; y (3) al hecho de que en muchas zonas se han aplicado medidas de restauración hidrológico forestal. Es el caso de lo acontecido en algunos sectores de los incendios de Las Peñas de Riglos (cod.: 31), Zuera y Castejón de Valdejasa (cod.: 34), tal y como se ha podido constatar a raíz de la consulta en la Dirección General de Gestión Forestal (Departamento de Agricultura, Ganadería y Medio Ambiente del Gobierno de Aragón). La heterogeneidad de las varianzas en cada intervalo, en especial en el I₂₀ que agrupa el mayor número de incendios, determinó que utilizáramos la prueba de Tamhane T2 para identificar la existencia de diferencias significativas entre intervalos.

Sin embargo, el gráfico de intervalos de confianza para la media que resume la distribución de los valores de fPAR dentro de las agrupaciones definidas, recogido en la Figura 2, pone de manifiesto la existencia de un patrón temporal en la distribución de los valores de fPAR a lo largo de los años de regeneración según la estación de los compuestos. En el caso de septiembre, se reconoce visualmente que no existen diferencias significativas entre los intervalos más alejados en el tiempo, aunque los valores son claramente superiores conforme aumenta la antigüedad, produciéndose un cambio de nivel entre los intervalos I₁₀-I₂₀ (entre 2 y 20 años de antigüedad) con los intervalos I₃₀-I₃₅ (< 30 años de antigüedad).

En este sentido, en incendios acontecidos en la década de los ochenta en el área prepirenaica occidental de Aragón, algunos de los cuales han sido considerados en este estudio, se ha descrito que, aunque 10-15 años no constituyen una franja temporal suficiente para que las formaciones vegetales afectadas se recuperen de manera integral, es suficiente tiempo para proporcionar buenos niveles de recubrimiento vegetal, así como cierta correspondencia fisionómico-florística entre las comunidades regeneradas y sus antecedentes (Pérez-Cabello, 2002). Indistintamente de la estrategia regenerativa de las comunidades afectadas, en estas primeras fases de regeneración son dominantes los estratos de vegetación herbácea y las leñosas inferiores a 0,5 m pudiendo cubrir por término medio el 50,3% de los nuevos espacios regenerados tras el fuego.

En el caso de mayo, los valores, como se ha dicho más arriba, presentan guarismos más elevados aunque se reproduce el mismo patrón temporal descrito para septiembre, excepto en el caso del I₂ que recoge valores inusuales y extraordinariamente altos. Ese último extremo se relaciona con la contribución de otras cubiertas (cultivos de cereal de secano) en los valores de LAI/fPAR, como consecuencia de la resolución espacial del producto MODIS y la estructura del paisaje de zonas quemadas situadas en la Depresión del Ebro (redes dendriformes de valles ocupados por cereal, rodeados de masas forestales). Por tanto, todo indica que, en zonas de gran heterogeneidad de cubiertas, los compuestos de otoño recogerían de manera más real el proceso de recuperación vegetal tras el fuego.

En la Figura 3 se representan los diagramas de dispersión correspondientes a las relaciones entre los valores de fPAR/LAI (compuestos de septiembre y mayo) y el tiempo trascurrido tras el fuego expresado en número de años. Los puntos representan a los valores medios y los errores típicos de las superficies quemadas cada año. Si bien los r² obtenidos son bajos debido a las peculiaridades de cada uno de los incendios considerados (entre 0,49 y 0,29), se identifica claramente una tendencia positiva en la que los valores de fPAR y, en menor medida LAI, aumentan en función del número de años tras el fuego. De tal modo que los incendios con mayor tiempo para su recuperación presentan en general valores más elevados. La correlación es más acusada en los compuestos de septiembre debido, especialmente, a su mejor ajuste en relación a los meses en los que suelen registrarse los incendios forestales, siendo muy similar la tendencia, independientemente de que se utilice fPAR o LAI.

Cabe preguntarse cómo se relacionan estos dos parámetros en las zonas quemadas. Físicamente la relación entre fPAR y LAI es exponencial inversa y muy dependiente del ángulo cenital solar (Jones & Vaughan 2010). En las zonas analizadas afectadas por incendio, los resultados, considerando cualquiera de los compuestos estacionales -mayo y septiembre- de las fechas disponibles y cuatro marcos geográficos, se identifica una relación de tipo lineal con r² superiores a 0,85 (Figura 4), debido probablemente al menor rango de los valores. Por tanto, aunque fPAR y LAI conceptualmente se refieren a variables biofísicas distintas del dosel vegetal, se encuentran estrechamente relacionadas (Zarate-Valdez et al., 2012; Glenn et al., 2008).

Si bien ambas variables son muy dependientes de las condiciones ambientales, LAI es una variable de mayor complejidad que se vincula de modo más directo a la capacidad fotosintética de las formaciones vegetales. Prueba de la mayor complejidad del LAI reside en el hecho de que mediante la inversión de la ecuación para la predicción de la radiación fotosintéticamente activa transmitida bajo el dosel (Norman & Campbell, 1989) es posible transformar fPAR en LAI siempre y cuando se disponga de datos de absortividad foliar y de los coeficientes de extinción de la radiación. Estos últimos dependen a su vez de las funciones de distribución angular de las hojas.

El análisis de la dinámica tras el fuego a mayor resolución temporal (año por año) se ha realizado para el periodo 2001-2010 en 6 incendios de la base de datos ocurridos en dicho intervalo temporal (Las Peñas de Riglos, CENAD (2), Zuera y Castejón de Valdejasa, Jaulín y Terrer), representativos además de los marcos geográficos contemplados en la zona de estudio. Para estos incendios se dispone de información sobre la situación previa y posterior a la perturbación, excepto en uno de los casos (Incendio de las Peñas de Riglos).

La Figura 5 y la Tabla 3 recogen la evolución temporal de los valores de LAI en el compuesto de septiembre. El análisis se ha restringido exclusivamente con este compuesto debido a su mayor proximidad a las fechas concretas del fuego, todos son incendios de verano, y al hecho de que el compuesto de mayo se ajusta peor al análisis de la evaluación de las repercusiones inmediatas del fuego en los valores de LAI. Del mismo modo, una vez analizada la correlación entre los valores de LAI y fPAR (apartado 3.1), el análisis se ha restringido sólo a los valores de LAI.

El primero de los gráficos contenidos en la Figura 5 (arriba a la izquierda) representa la trayectoria multianual del incendio de las Peñas de Riglos de 2001 (código 31), perteneciente al marco geográfico de las Sierras occidentales del Prepirineo oscense, dominado en este caso por comunidades de carácter submediterráno (Quercus gr. cerrioides y Buxus sempervirens), plantaciones de Pinus nigra y matorrales.

En este caso, la curva sólo recoge el proceso de recuperación de los valores de LAI transcurridos dos años tras el fuego debido a que el incendio tuvo lugar en 2001 y la serie de datos que se maneja comienza en 2003, ya que el compuesto utilizado combina la información espectral procedente de dos satélites, Terra y Aqua, siendo este último lanzado más tarde, en 2002. En general la curva muestra gran estabilidad interanual, a pesar de la posible influencia en los valores de los trabajos de apeo y extracción de la madera quemada y de las tareas de restauración de la vegetación, a las que se vio sometido desde el año siguiente al fuego; al menos en los 2/3 de la superficie afectada (2.212 ha). No obstante, destaca el descenso de los valores en 2005, como consecuencia de los efectos de la sequía del año hidrológico 2004-2005. A partir de un simple análisis de la serie mensual de datos de precipitación de la estación de Huesca (Red de Estaciones Agro-meteorológicas del SIAR –MAGRAMA-), próxima al incendio de las Peñas de Riglos, se recoge en ese año una reducción del 33% respecto a la media histórica.

El gráfico con el código 32 representa la evolución del incendio de 2004 que tuvo lugar dentro del CENAD, dominado por formaciones de carácter subarbustivo mediterráneo (Rosmarinus officinalis, Genista scorpius, etc.) con pastizales xerofíticos. En la curva se advierte el descenso provocado por el fuego (decremento del 57% del valor previo de LAI) y la asténica recuperación durante los años siguientes cuyo inicio hay que situarlo en 2006, dos años después del fuego. Dicho retraso podría responder a las situaciones de estrés hídrico de 2005 en esta zona. Este mismo proceso de ralentización de la regeneración con motivo de la menor cantidad de precipitación en 2005 fue descrito en la estación experimental de Peñaflor (EEP) (cercana CENAD), en donde se realizó un fuego controlado sobre formaciones vegetales similares a las de la zona del CENAD (Pérez-Cabello et al., 2012), recogiéndose un 25% menos de precipitación media anual.

El gráfico con el código 34 se corresponde con el incendio de Zuera y Castejón de Valdejasa de 2008, que afectó en su mayoría a una masa de pinar de Pinus halepensis con Q. coccifera. El valor de 2008 recoge el impacto del fuego, mientras que los registros de 2009 y 2010 representan la recuperación de los valores de LAI. El incremento de los valores de LAI en este caso puede explicarse por la eficacia de la estrategias activas de respuesta al fuego que detentan las dos especies principales afectadas: dispersión de propágulos en el caso de Q. coccifera, y germinación en el caso de P. halepensis, ya documentadas ampliamente en zonas incendiadas próximas al incendio de 2008 (Vicente-Serrano et al., 2011; Tanase et al., 2011).

Por otro lado, en este incendio hay que destacar la alta variabilidad que recogen los valores de LAI en los años previos al fuego, así como el inexplicable descenso de 2003. En este sentido, la singular estructura del paisaje de este sector, caracterizada por sistemas dendriformes de estrechos valles ocupados por cultivos de cereal de invierno e interfluvios cubiertos por arboledas de P. halepensis, podría explicar dicha variabilidad, en relación también con la resolución espacial de los productos MODIS. A pesar de que en septiembre la campaña cerealista ha finalizado, el diferente estado de los cultivos en relación a los restos de la siega y la colonización de vegetación natural, podría influir de manera significativa.

La fila inferior de la figura 5 recoge la evolución en tres incendios que tuvieron lugar en 2009. El impacto del fuego queda registrado en el acusado descenso en los valores de LAI_sep en el compuesto de 2009. Antes del fuego, destaca la estabilidad de los valores de LAI_sep, sólo interrumpida por un leve repunte en el año 2004. En los tres casos, a pesar de partir de niveles previos al fuego diferentes, como consecuencia de que representan marcos geográficos distintos (matorrales suarbustivos mediterráneos y pastizales xerofíticos -38-, pinares de halepensis -36- y matorrales acidófilos -37), los valores tras el fuego presentan magnitudes similares, en torno a 0,20, siendo significativamente diferentes a los de los años previos al fuego, mostrando un descenso del orden del 75%.

Por tanto, las líneas de evolución presentadas a esta escala temporal reproducen los niveles de impacto y recuperación del índice de área foliar de la vegetación tras el fuego en diferentes formaciones vegetales. El descenso de los valores tras el fuego puede aportar nuevas claves de interpretación biofísica en relación con las mediciones de la severidad, aspecto muy importante desde el punto de vista de la recuperación de las zonas quemadas. Asimismo, las magnitudes en el incremento permiten cuantificar y determinar espacialmente el grado de recomposición de las funciones eco-fisiológicas de la vegetación en su relación con la atmósfera.

Por otro lado, un análisis más exhaustivo de diferentes factores que hayan contribuido de modo significativo en los valores de LAI/fPAR a lo largo del tiempo (termo-pluviometría, severidad del fuego, estrategia reproductiva), así como de información de campo, proporcionaría nuevos elementos de análisis para interpretar, de modo más preciso, los valores y los patrones espacio/temporales de la información derivada del producto MCD13A2 de MODIS.