Validacion-de-un-modelo-de-perturbacion-forestal-por-uso-de-lena-en-Mexico

Vista previa del material en texto

UNIVERSIDAD NACIONAL AUTÓNOMA DE MÉXICO
PROGRAMA DE MAESTRÍA Y DOCTORADO EN GEOGRAFÍA
VALIDACIÓN DE UN MODELO DE PERTURBACIÓN FORESTAL
POR USO DE LEÑA EN MÉXICO
TESIS
QUE PARA OPTAR POR EL GRADO DE:
MAESTRO EN GEOGRAFÍA
PRESENTA:
ALEXANDER QUEVEDO CHACÓN
DIRECTORES DE TESIS:
Dr. ADRIÁN GHILARDI (CIGA - UNAM)
Dra. YAN GAO (CIGA - UNAM)
Ciudad Universitaria, Cd. Mx. Enero 2018
 
UNAM – Dirección General de Bibliotecas 
Tesis Digitales 
Restricciones de uso 
 
DERECHOS RESERVADOS © 
PROHIBIDA SU REPRODUCCIÓN TOTAL O PARCIAL 
 
Todo el material contenido en esta tesis esta protegido por la Ley Federal 
del Derecho de Autor (LFDA) de los Estados Unidos Mexicanos (México). 
El uso de imágenes, fragmentos de videos, y demás material que sea 
objeto de protección de los derechos de autor, será exclusivamente para 
fines educativos e informativos y deberá citar la fuente donde la obtuvo 
mencionando el autor o autores. Cualquier uso distinto como el lucro, 
reproducción, edición o modificación, será perseguido y sancionado por el 
respectivo titular de los Derechos de Autor. 
 
 
 
Resumen
La validación de modelos espacialmente expĺıcitos está tomando cada vez mayor relevancia,
ante el creciente número de propuesta de modelación de los fenómenos ambientales, y la nece-
sidad de los usuarios de conocer su grado de confiabilidad. La presente tesis busca determinar
la confianza en un modelo espacio temporal que simula los efectos de la recolección de leña
sobre la vegetación (MoFuSS). Para este fin se empleó un análisis de series de tiempo de ı́ndi-
ces de vegetación obtenidas de MODIS para evaluar la coincidencia entre áreas que ha sufrido
perturbaciones forestales y las predicciones del modelo MoFuSS, para los estados de Yucatán,
Michoacán y Jalisco en el periodo comprendido entre 2008 y 2016. Los resultados obtenidos
indican que para los dos primeros estados existe una coincidencia espacial entre las perturba-
ciones detectadas mediante MODIS y las estimaciones realizadas por el modelo MoFuSS. En
el caso de Jalisco, no se encontró coincidencia espacial. La coincidencia espacial mostrada en
Yucatán y Michoacán indican que es posible tener confianza en las predicciones hechas por el
modelo MoFuSS en cuanto a las áreas que presentan perturbación, siendo que éstas pueden ser
ocasionadas por diversos factores, por lo que no se puede atribuir dicha perturbación a un solo
factor. El análisis de series de tiempo de imágenes satelitales mostró ser de gran utilidad en el
proceso de validación de este modelo. Se recomienda complementar esta técnica con estudios
de campo y trabajar a un resolución más detallada, usando sensores como Landsat 8 y Sentinel
2.
3
Agradecimientos
Institucionales
Al Consejo Nacional de Ciencia y Tecnoloǵıa (CONACYT) por la beca otorgada durante la
realización de la maestŕıa.
Al Centro de Investigaciones en Geograf́ıa Ambiental (CIGA) de la Universidad Nacional
Autónoma de México (UNAM).
A mis asesores el Dr. Adrián Ghilardi y a la Dra. Yan Gao, por su constante apoyo y
apertura durante el desarrollo de este proyecto. Al Dr. Jean Francois Mas Caussel por su ayuda
y sus acertadas observaciones.
Al comite sinodal, Dr. Ernesto Vega, Dr. Alberto Gomez-Tagle y al Dr. Robert Bailis.
Investigación realizada gracias al Programa UNAM-DGAPA-PAPIIT No IA101513: “Análi-
sis geoespacial de la degradación forestal por la extracción de madera para leña y carbón vegetal
en el centro de México.” Ene 2013 – Dic 2014.
Investigación realizada gracias al Fondo Sectorial CONACYT-SENER Sustentabilidad energéti-
ca Proyecto No 246911. Clúster de Biocombustibles Sólidos para la Generación Térmica y
Eléctrica (BCS-GTE). Nov 2016 – Oct 2020.
Investigación realizada gracias al Programa de Apoyo de Proyectos de Investigación e Inova-
ción Tecnológica PAPIIT IA104117, UNAM. Mapeo de la degradación forestal en México uti-
lizando datos de serie de tiempo de ı́ndice de vegetación del sensor MODIS.
4
Personales
A mis padres y mi hermana por todo su apoyo en la distancia. A Adri por todo su amor,
ternura y comprensión.
Ale Larrazabal, Renato y Paty por acogerme desde mi llegada a México.
A Adrián Ghilrdi por su sinceridad y apoyo mucho más allá de lo académico.
A nuestro buen amigo Domingo por toda su generosidad; al sector 7G especialmente a Jose
Luis y al Mike Salinas; a Abdul Reyes (Princess) por toda su solidaridad; a mis compañeros de
maestŕıa y a todos los que en la cercańıa o la distancia han estado.
Y por supuesto a la vida que me ha dado y enseñado tanto.
5
Índice general
1. Introducción 12
1.1. Pregunta de investigación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
1.2. Objetivo General . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
1.3. Objetivos espećıficos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
2. Marco Teórico 15
2.1. Bosque y perturbaciones forestales . . . . . . . . . . . . . . . . . . . . . . . . . . 15
2.2. Dendroenerǵıa . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
2.3. Modelos ambientales . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.4. Monitoreo de la cubierta vegetal . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
3. Metodoloǵıa 22
3.1. Área de estudio . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
3.2. Datos y Sofware . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
3.3. Imagenes MODIS - NDVI . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3.4. Preprocesamieto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
3.4.1. Filtros de calidad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
3.4.2. Interpolación y Filtro Savitzky Golay . . . . . . . . . . . . . . . . . . . . 29
3.5. Bfast monitor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
3.5.1. Modelo Tendencial y Estacional . . . . . . . . . . . . . . . . . . . . . . . . 32
3.5.2. Monitoreo del cambio . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
6
3.5.3. Peŕıodo de histórico estable . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.5.4. Bfast Spatial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
3.6. Modelización de escenarios de ahorro de leña - MoFuSS . . . . . . . . . . . . . . 35
3.6.1. Mapas de fricción . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
3.6.2. Componente IDW . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
3.6.3. Demanda y Oferta de leña . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
3.6.4. Pérdida y ganancia de bosque. . . . . . . . . . . . . . . . . . . . . . . . . 39
3.7. Análisis de puntos calientes y puntos fŕıos . . . . . . . . . . . . . . . . . . . . . . 40
3.8. Validación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.8.1. Análisis de correlación . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.8.2. Contraste de hipótesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.8.3. Índice de coincidencia difusa . . . . . . . . . . . . . . . . . . . . . . . . . 42
4. Resultados y discusión 43
4.1. Estimación de pérdida de la cobertura vegetal . . . . . . . . . . . . . . . . . . . . 43
4.2. Impactos proyectados por extracción de leña a partir de MoFuSS . . . . . . . . . 43
4.3. Agrupaciones de áreas con pérdida de NDVI (Puntos calientes - puntos fŕıos) . . 47
4.4. Resultados de las estrategias de validación de MoFuSS . . . . . . . . . . . . . . . 49
4.4.1. Análisis de correlación entre Biomasa No Renovable (NRB) y pérdida del
NDVI . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
4.4.2. Comparación estad́ıstica . . . . . . . . . . . .. . . . . . . . . . . . . . . . 49
4.4.3. Índice de coincidencia difusa entre áreas de NRB y Agrupaciones de áreas
con pérdida de NDVI . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
4.5. Discusión . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
4.5.1. Validación de MoFuSS usando series de tiempo de NDVI . . . . . . . . . 59
4.6. Conclusiones . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
4.7. Recomendaciones . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
7
Índice de figuras
2-1. Ciclo de desarrollo de un modelo . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
3-1. Flujo de trabajo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
3-2. Área de estudio . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3-3. Algoritmo de filtrado empleando la capa de calidad de NDVI. . . . . . . . . . . . 29
3-4. Representación gráfica del periodo histórico y de monitoreo empleado para el
Análisis de Serie de Tiempo con Bfast. . . . . . . . . . . . . . . . . . . . . . . . . 31
3-5. Ejemplo de la salida gráfica del algoritmo Bfast monitor aplicado en una zona
de selva baja subcaducifolia en Michoacán. . . . . . . . . . . . . . . . . . . . . . 35
3-6. Algoritmo MoFuSS simplificado . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
4-1. Pérdida sobre la cubierta vegetal en unidades de NDVI para el estado de Yucatán
(2008 - 2016) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
4-2. Pérdida sobre la cubierta vegetal en unidades de NDVI para el estado de Mi-
choacán (2008 - 2016) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
4-3. Pérdida sobre la cubierta vegetal en unidades de NDVI para el estado de Jalisco
(2008 - 2016) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
4-4. Estimacion de la extracción de biomasa leñosa a tasas que exceden la regeneración
natural (NRB) dentro del peŕıodo 2008 - 2016, para el estado de Yucatán. . . . . 45
4-5. Estimación de la extracción de biomasa leñosa a tasas que exceden la regeneración
natural (NRB) dentro del peŕıodo 2008 - 2016, para el estado de Michoacán. . . . 46
8
4-6. Estimación de la extracción de biomasa leñosa a tasas que exceden la regeneración
natural (NRB) dentro del peŕıodo 2008 - 2016, para el estado de Jalisco. . . . . . 46
4-7. Puntos fŕıos y puntos calientes de pérdida de NDVI, calculados mediante el ı́ndice
Getis-Ord para el estado de Yucatán . . . . . . . . . . . . . . . . . . . . . . . . . 47
4-8. Puntos fŕıos y puntos calientes de pérdida de NDVI, calculados mediante el ı́ndice
Getis-Ord para el estado de Michoacán . . . . . . . . . . . . . . . . . . . . . . . . 48
4-9. Puntos fŕıos y puntos calientes de pérdida de NDVI, calculados mediante el ı́ndice
Getis-Ord para el estado de Jalisco . . . . . . . . . . . . . . . . . . . . . . . . . . 48
4-10. Correlaciones locales entre los valores de NRB y la perdida de NDVI . . . . . . . 50
4-11. Diagramas de cajas para cada uno de los grupos establecidos . . . . . . . . . . . 52
4-12. Comparación NRB generado por MofuSS (panel izquierdo) y NRB aleatorio (pa-
nel derecho) para el estado de Yucatán . . . . . . . . . . . . . . . . . . . . . . . . 54
4-13. Comparación NRB generado por MofuSS (panel izquierdo) y NRB aleatorio (pa-
nel derecho) para el estado de Michoacán . . . . . . . . . . . . . . . . . . . . . . 54
4-14. Comparación NRB generado por MofuSS (panel izquierdo) y NRB aleatorio (pa-
nel derecho) para el estado de Jalisco . . . . . . . . . . . . . . . . . . . . . . . . . 55
4-15. Esquema de la comparación espacial entre las áreas de NRB y NRB aleatorio
con las áreas definidas como clusters de perdida. . . . . . . . . . . . . . . . . . . 55
4-16. Relación entre el ı́ndice coincidencia difusa (eje Y) y las comparaciones entre el
NRB de MoFuSS y los clústers de pérdida de NDVI, y entre el NRB aleatorio
(NRBa) y los clústers de pérdida de NDVI; en función de la distancia (eje X),
para el estado de Yucatán. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
4-17. Relación entre el ı́ndice coincidencia difusa (eje Y) y las comparaciones entre el
NRB de MoFuSS y los clústers de pérdida de NDVI, y entre el NRB aleatorio
(NRBa) y los clústers de pérdida de NDVI; en función de la distancia (eje X),
para el estado de Michoácan. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
9
4-18. ÍRelación entre el ı́ndice coincidencia difusa (eje Y) y las comparaciones entre el
NRB de MoFuSS y los clústers de pérdida de NDVI, y entre el NRB aleatorio
(NRBa) y los clústers de pérdida de NDVI; en función de la distancia (eje X),
para el estado de Jalisco. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58
4-19. Semivariograma de NRB ajustado a un modelo esférico para el estado de Yucatán118
4-20. Semivariograma de NRB ajustado a un modelo esférico para el estado de Michoacán119
4-21. Semivariograma de NRB ajustado a un modelo esférico para el estado de Jalisco 119
10
Índice de cuadros
3-1. Datos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
3-2. Descripción de la capa de calidad de los ı́ndices de vegetación. . . . . . . . . . . . 27
3-2. Descripción de la capa de calidad de los ı́ndices de vegetación. . . . . . . . . . . . 28
3-2. Descripción de la capa de calidad de los ı́ndices de vegetación. . . . . . . . . . . . 29
3-3. Descripcion de las fuentes de oferta de leña. . . . . . . . . . . . . . . . . . . . . . 38
4-1. Correlaciones globales entre los valores de NRB y la perdida de NDVI. . . . . . . 49
4-2. Autocorrelación espacial. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
4-3. Prueba de normalidad a los valores de NRB. . . . . . . . . . . . . . . . . . . . . 51
4-4. Estad́ısticas descriptivas para cada uno de los grupos establecidos. . . . . . . . . 51
4-5. Comparación entre los tres grupos definidos empleando la prueba de Kruskal-
Wallis. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
4-6. Comparación entre los valores de NRB al interior de los clusters y los valores
fuera de los cluesters, empleando la prueba U de Mann-Whitney. . . . . . . . . . 53
11
Caṕıtulo 1
Introducción
El uso tradicional de leña y carbón vegetal en los páıses en desarrollo es a menudo identificado
como una causa de degradación forestal, e inclusive de deforestación [Hosonuma et al., 2012]; lo
que a su vez tiene implicaciones ambientales por la pérdida de bienes y servicios ecosistémicos
y por las emisiones netas de gases de efecto invernadero [Bailis et al., 2015]. De igual manera,
la escasez de estos recursos pone en riesgo la seguridad energética a un tercio de la población
mundial, que depende exclusivamente de la leña para cocinar sus alimentos y calentar sus
hogares [Legros et al., 2009]. Se espera que la demanda global de leña y carbón vegetal continúe
en aumento (en valores absolutos) al menos hasta el año 2030 [Bonjour et al., 2013; Riahi et al.,
2012]. Bajo prácticas de manejo sustentable, la leña y el carbón vegetal representan una fuente
de enerǵıa renovable y de ingresos para la población rural y urbana de bajos recursos [Masera
et al., 2015].
En el ámbito de las poĺıticas públicas se necesitan evaluaciones robustas y precisas que
tengan en cuenta los efectos espaciotemporales de la cosecha de leña para mejorar el impacto
de las intervenciones, tales como programas de cocinas y hornos de carbón mejorados. En el
pasado los aciertos de estos programas se han asumido como una cuestión de confianza en la
tecnoloǵıa y no a partir de evaluaciones cient́ıficas. Por lo tanto existe una fuerte necesidad
de modelos que permitan evaluar de manerarobusta las decisiones de manejo de este recurso
[Ghilardi et al., 2016].
12
En este contexto entre los años 2011 y 2015, Ghilardi y colaboradores (2016) desarrollaron
el modelo MoFuSS (Modeling fuelwood savings scenario, por sus siglas en inglés), el cual es
espacialmente expĺıcito y dinámico, y simula el efecto de la recolección de leña en la vegetación
local. Este modelo se desarrolló siguiendo tres preguntas gúıas: (1)¿Cuánta leña se cosecha en un
lugar determinado dentro de un marco de tiempo espećıfico?, (2) ¿Cómo responde la vegetación
a esta presión, medida con la existencia de biomasa y las tasas de crecimiento?, y (3)¿Cómo son
los cambios en la demanda de leña (por ejemplo, a través de la implementación de estufas que
ahorran combustible), en términos de la alteración en los patrones de cosecha y recrecimiento
con el tiempo? [Ghilardi et al., 2016]. Estas preguntas definen el dominio del modelo y por ende
la aplicabilidad del mismo.
MoFuSS fue diseñado para ser aplicado en áreas donde la leña es una fuente de enerǵıa
importante para el sector residencial, en donde este recurso es recolectado para uso propio y
para venta local, y no para áreas en que los patrones de cosecha cumplen fines comerciales de
gran escala [Masera et al., 2015]. Este modelo fue aplicado en 2015 en cuatro departamentos de
Honduras, donde encontraron, mediante juicios de expertos, congruencia entre los resultados
del modelo y los patrones esperados basados en el conocimiento de 30 años de estudios, sobre
la recolección y cosecha de leña con fines residenciales [Ghilardi et al., 2016].
Dentro los módulos que componen el modelo MoFuSS se encuentra uno destinado a integrar
la incertidumbre inherente de algunos parámetros de entrada (i.e. crecimiento de biomasa, pro-
porción de zonas visitadas al menos una vez para la cosecha de leña). Sin embargo este modelo
no ha sido validado empleando datos independientes [Ghilardi et al., 2016]. La importancia
de la validación de los modelos radica en que ésta conlleva una mejora de la comprensión de
la capacidad del modelo para simular con precisión los cambios [van Vliet et al., 2016]. La
validación de los modelos ambientales permite que éstos se sumen a los avances cient́ıficos y
a la aceptación por parte de los actores sociales involucrados y los tomadores de decisiones,
quienes se encuentran cada vez más interesados en comprender la incertidumbre de estos mo-
delos [O’Hagan, 2012; van Vliet et al., 2016]. Es por esto que esta investigación se centra en la
evaluación de MoFuSS, empleando datos independientes obtenidos a partir de series de tiempo
13
de imágenes MODIS, tomando como caso de estudio tres estados del territorio mexicano.
1.1. Pregunta de investigación
¿Existe coincidencia espacial entre la distribución de las zonas susceptibles a la perturbación
forestal reportadas en el modelo, y los cambios observados sobre la cubierta forestal a partir de
los ı́ndices de vegetación MODIS de 2008 a 2016?
1.2. Objetivo General
Validar la distribución espacial de zonas con diferente susceptibilidad a la perturbación
forestal por uso de leña, reportadas en el modelo “Modeling Fuelwood Savings Scenario” (Mo-
FuSS).
1.3. Objetivos espećıficos
Estimar tendencias espacio temporales de pérdida de cobertura vegetal en México utili-
zando ı́ndices de vegetación de MODIS entre el periodo 2000-2015; y reclasificarlas en un
ı́ndice ordinal.
Estimar la degradación forestal por extracción y uso de leña, mediante el modelo “Mode-
ling fuelwood savings scenario”.
Comparar los resultados obtenidos de los análisis anteriores para evaluar la correspon-
dencia espacial.
14
Caṕıtulo 2
Marco Teórico
2.1. Bosque y perturbaciones forestales
Existen muchas definiciones de lo que se considera un bosque como exponen Chazdon
et al. [2016], algunas de estas definiciones son adoptadas por instituciones gubernamentales
y académicas como marco de referencia para el estudio y generación estrategias de manejo y
conservación de los bosques. Tal es el caso de la definición empleada por la Organización de las
Naciones Unidas para la Alimentación y la Agricultura-FA0; de forma inicial este trabajo se
ajusta a esta definición:
”Tierras que se extienden por más de 0,5 hectáreas dotadas de árboles de una altura
superior a 5 metros y una cubierta de dosel superior al 10 por ciento, o de árboles
capaces de alcanzar esta altura in situ. No incluye la tierra sometida a un uso
predominantemente agŕıcola o urbano”[FAO, 2012, pág.3]
Dadas las problemáticas ambientales actuales [IPCC, 2014] y los cambios en la cubierta te-
rrestre que afectan a los bosques en todo el mundo, se ha hecho crucial comprender las dinámicas
de los bosques para su conservación [Anderson-Teixeira et al., 2015]. Parte de estas dinámicas
son las perturbaciones antrópicas tales como la deforestación y la degradación forestal, las cuales
implican pérdidas de biomasa.
En este orden de ideas el presente trabajo emplea la definición de deforestación dada por
15
la FAO [2012, pág.6]: ”la conversión de los bosques a otro tipo de uso de la tierra o la reduc-
ción permanente de la cubierta de dosel, por debajo del umbral mı́nimo del 10 por ciento”. Es
importante aclarar que esta definición excluye las áreas en las cuales los arboles son extráıdos
debido al aprovechamiento forestal y en las cuales se espera que el bosque se regenere de manera
natural o bajo manejo silv́ıcola.
La degradación forestal por otra parte cuenta con un numero de definiciones más amplio,
y que en ocasiones son dif́ıciles de hacer operativas para su evaluación [Ghazoul et al., 2015;
Simula, 2011]. Tal es el caso de la definición empleada por la FAO [2012, pág.28]: ”La degradación
forestal es la reducción de la capacidad del bosque de proporcionar bienes y servicios.”, por esta
razón se decidió tomar como marco de referencia la definición aportada por Ghazoul et al. [2015,
pág.630]:
”La degradación es un estado de sucesión detenida o arrestada, que se establece en
contraste con los sistemas naturales, los cuales experimentan transiciones constantes
sujetas tanto a los efectos de vecindad como a las perturbaciones naturales. El
corolario de esto es que un bosque no se degrada siempre y cuando mantenga una
dinámica que facilite la recuperación a uno de los estados estables anteriores. Una
vez que se presenta un estado de sucesión detenida, se requiere intervención externa
para recuperar trayectorias sucesionales, para lo cual pueden estar vinculados los
procesos de toma de decisiones y de gestión.”
Dicha definición presenta ventajas respecto a la definición de la FAO, ya que abre mayores
posibilidades para la estimación de la degradación mediante técnicas de teledetección y por
ende amplia las posibilidades de manejo.
2.2. Dendroenerǵıa
Para el año 2011, el 80% de la maderera extráıda en todo el mundo se usaba con fines
energéticos [FAO, 2011], esto implica que más de 2,700 millones de personas dependen de esta
fuente de enerǵıa para cocinar y calentarse, y se espera que para el año 2030 el numero de
16
usuarios se incrementará a 2,800 millones de personas. Los hogares consumidores de leña se
localizan principalmente en los páıses en v́ıa de desarrollo [FAO, 2011; IEA, 2010].
La extracción de leña puede contribuir a la degradación de los bosques o la deforestación, si
la madera es extráıda más rápido de lo que el bosque puede regenerarse (biomasa no renovable).
Estos procesos de cambio impactan negativamente sobre las condiciones ambientales a escala
local y generan emisiones de CO2 [Bailis et al., 2015; Masera et al., 2015]. Se estima que el mal
manejo de este recurso (cosecha insostenible y combustión incompleta), contribuye a cerca del
2% de las emisiones totales de gases de efecto invernadero a nivel mundial.
A pesar de la gran importancia de esterecurso se carece de estad́ısticas detalladas de
consumo y proyecciones futuras, en gran medida debido a que una alta proporción de la cosecha
de leña se realiza en zonas rurales para uso propio y mercados locales. Adicionalmente en
algunos páıses la cosecha de leña se realiza de manera ilegal, haciendo dif́ıcil la adquisición
de estad́ısticas fiables. Un tercer motivo es que las estad́ısticas generalmente se encuentran
dispersas entre diferentes entidades [UNECE/FAO, 2012].
En el año 2016 se desarrolló Mofuss, ante la necesidad de obtener una cuantificación de la
biomasa no renovable a través de un método robusto y ampliamente replicable, que permita
orientar proyectos de compensación de carbono para cocinas, inventarios nacionales y estrategias
sostenibles de manejo forestal [Ghilardi et al., 2016].
2.3. Modelos ambientales
La modificación humana al ambiente es un proceso complejo y multidimensional cuya com-
prensión requiere el conocimiento de una amplia gama de disciplinas cient́ıficas. Este proceso
es el resultado directo de la toma de decisiones por parte de los seres humanos, teniendo una
gran diversidad de causas, que implican varios factores, aśı como entornos/escalas regionales y
globales [Magliocca et al., 2015].
En los estudios ambientales, las principales razones para la modelación son la generación
y el intercambio de conocimiento para apoyar la toma de decisiones con fines de gestión o de
17
Figura 2-1: Ciclo de desarrollo de un modelo, en cual se indican las etapas y los procesos de
validación (Magliocca 2015). El ciclo inicia con la identificación de la entidad del problema, a
partir del cual se desarrollan el modelo conceptual y el computacional y finaliza con la aplicación
del modelo.
desarrollo de poĺıticas publicas estratégicas. El modelado es un proceso cient́ıfico muy común
de simplificación de la realidad para mejorar la comprensión del sistema, por lo cual en las
ultimas dos décadas se han desarrollado y aplicado una amplia gama de modelos, con estos
fines [Jakeman et al., 2006, 2008; van Vliet et al., 2016].
El desarrollo del modelo de un sistema se puede describir como una proceso secuencial de
varios pasos, no necesariamente lineal, sino mas bien iterativo, donde las iteraciones resultan
esenciales para su generación [Jakeman et al., 2006; Sargent, 2005; Magliocca et al., 2015; van
Vliet et al., 2016]. A la luz de esto Magliocca et al. [2015] propone un esquema simplificado del
ciclo de vida de un modelo (Figura.2-1), a partir del esquema de Robinson et al. [2006].
A continuación se describen brevemente los pasos que componen el ciclo de desarrollo del
18
modelo, haciendo énfasis en la validación operacional del mismo, debido a que el objetivo de
esta investigación se centra en este paso.
El punto de inicio en el desarrollo de un modelo según Magliocca et al. [2015], Figura 2-1,
es la entidad del problema, se puede definir como el fenómeno que es objeto de investigación.
Este paso incluye la selección de posibles variables explicativas, procesos involucrados y los
ĺımites del sistema a modelar basándose en preguntas de investigación [Sargent, 2005; van Vliet
et al., 2016]. Durante el desarrollo del modelo y con el aumento en la compresión del sistema
es altamente probable que el modelador regrese a este punto [Jakeman et al., 2006].
El conocimiento de la identidad del problema permite dar el siguiente paso en el ciclo pro-
duciendo el modelo conceptual. En este punto se describe la problemática de investigación en
términos de componentes y relaciones (e.g. Factores de control o determinantes de deforestación
o de cambio de uso del suelo). La finalidad del modelo conceptual es hacer expĺıcito como se en-
tiende el sistema, en la práctica esto puede expresarse mediante una representación matemática,
lógica o gráfica [Sargent, 2005; van Vliet et al., 2016].
La validación conceptual evalúa si el modelo conceptual y los procesos, conceptos y supuestos
que en éste se representan, son adecuados y lógicos en el contexto del propósito del modelo [van
Vliet et al., 2016].
El modelo computacional es la traducción del modelo conceptual en el código informático.
Esta implementación puede ser verificada, ya que es posible comprobar exactamente cuán bien
el código de software representa cada componente del modelo conceptual [Sargent, 2005].
Según el ciclo presentado Magliocca et al. [2015] la aplicación del modelo se deifne como una
instancia del modelo computacional que incluye los datos y parámetros para un caso particular.
El proceso para ajustar los parámetros y mejorar el desempeño de un modelo se denomina
calibración [van Vliet et al., 2016].
La validación operativa del modelo se define como la determinación de que los resultados de
salida del modelo tengan un rango satisfactorio para el propósito previsto[Sargent, 2005]. Los
modelos ambientales rara vez pueden ser evaluados completamente, principalmente debido a la
heterogeneidad de sus datos y a la gama de factores que influyen en la utilidad de sus productos
19
[Jakeman et al., 2006]. Adicionalmente la validación es un proceso costoso y requiere mucho
tiempo para determinar que un modelo es absolutamente válido sobre el dominio completo de
su aplicabilidad predeterminada. Debido a esto, en su lugar se realizan pruebas y evaluaciones
hasta que se tiene la suficiente confianza, y el modelo puede considerarse valido para aplicación
para la que fue previsto [Sargent, 2005].
En cuanto a los métodos de evaluación o validación de modelos ambientales no existe un
acuerdo universal sobre un umbral en la bondad de ajuste para considerar un modelo válido o
no [Pontius et al., 2004; Jakeman et al., 2006; Uusitalo et al., 2015], sino que este criterio debe
ser definido en función del propósito del modelo.
Referente a las técnicas de validación se recomienda no emplear aquellas que brinden un
solo número para expresar el ajuste, ya que éstas dif́ıcilmente proporcionarán información sobre
cómo mejorar el modelo [Pontius et al., 2004]. En el componente espacial la escala es de crucial
importancia en la evaluación debido a que los resultados pueden ser sensibles a ciertos patrones,
que solo son evidentes en ciertas escalas. [Pontius et al., 2004; Uusitalo et al., 2015].
2.4. Monitoreo de la cubierta vegetal
Las técnicas de percepción remota son una herramienta exitosa para detectar el cambio sobre
la superficie terrestre y en la actualidad juegan un papel clave en la estimación de cambios sobre
la cubierta vegetal. Esto se debe a dos factores fundamentales: la creciente disponibilidad de
repositorios de imágenes satélitales con diferentes resoluciones espaciales y temporales, y el
desarrollo de diversos procedimientos de monitoreo para detectar los cambios [Mayaux et al.,
2008; Mas, 1999; Kuenzer et al., 2015].
Uno de los procedimientos que está dando mayores resultados en la detección de cambios
es el análisis de series de tiempo; éste emplea métodos estad́ısticos para analizar y modelar
observaciones equidistantes en el tiempo [Kuenzer et al., 2015]. El análisis de series de tiempo
para la percepción remota se entiende como el monitoreo temporalmente denso de la dinámica
de la superficie terrestre en un periodo definido.
20
Este método presenta ventajas sobre el método tradicional, en el cual se evalúan dos fechas
ya que permite conocer con mayor precisión temporal la ocurrencia del cambio; adicionalmente
permite establecer si el cambio se produce debido a una variabilidad en la estacionalidad o es
producto de un cambio abrupto en la tendencia de la cubierta[Verbesselt et al., 2010a].
Los datos procedentes de los sensores pueden ser series de valores de reflectancia o valores
derivaros como el Índice Normalizado de Vegetación (NDVI, por sus siglas en ingles). La idea
básica detrás de los diferentes ı́ndices de vegetación existenteses que la combinación algebraica
de bandas espectrales detectadas de forma remota, tiene la capacidad de revelar información
de gran utilidad sobre el estado de cubierta vegetal, la capacidad fotosintética, densidad y
distribución foliar, contenido de agua en las hojas, deficiencias minerales, entre otros datos
relevantes para diferentes estudios de monitoreo de la cubierta vegetal. [Yengoh et al., 2015b].
21
Caṕıtulo 3
Metodoloǵıa
Este capitulo describe la metodoloǵıa de esta investigación. La Figura 3-1 representa de
manera simplificada el flujo de trabajo que se siguió para alcanzar el objetivo planteado. Esta
metodoloǵıa es la fusión de métodos y técnicas empleadas en estudios previos de monitoreo de
la vegetación y de modelos de cambio de cobertura y uso del suelo.
3.1. Área de estudio
El área de estudio se seleccionó con base en el trabajo realizado por Ghilardi [2008], en el cual
se desarrolla un ı́ndice de prioridad por consumo de leña. En la presente tesis se seleccionaron
estados en cuyos municipios predomina una de las categoŕıas este ı́ndice (alto, medio y bajo),
siendo Yucatán, Michoacán de Ocampo y Jalisco los estados seleccionados (Figura 3-2), los
cuales se describen brevemente a continuación:
Yucatán: Estado localizado en el sur de México, con un área de 39,524 km2. Según el
Instituto Nacional de Estad́ıstica y Geograf́ıa INEGI [2013] la cobertura predominante
en este estado es la vegetación de selva mediana caducifolia con un 41% del área total,
seguida de pastizales inducidos y cultivados con 18%; en tercer lugar se encuentra la
vegetación de selva mediana caducifolia con 15%, el restante 26% es ocupado por otras
25 categoŕıas. En cuanto a población, en el año 2010 [INEGI, 2010] Yucatán contaba con
22
Figura 3-1: Flujo de trabajo, en el cual se detallan los datos, procesos y resultados de la
metodologia utilizada.
23
/ 
Imagenes / 
MODIS 
Construcción de 
mosaicos, cambio 
de sistema de 
referencia 
Aplicación de un 
filtro de calidad 
Aplicación un filtro 
de señales 
Detección iterativa 
de cambio sobre 
la cobertura 
forestal 
/
Area de perdida~~ 
sobre la cobertura 
forestal 
Estadistica 
Getis-Ord Gi* 
,"cluesters " de 
cambio 
!rametros de enttrada ie Mofuss, caso México 
I Mofuss I 
/
Biomasa no Renovable / 
(NRB) 
~ Análisis de correlación 11--------1 
~--If....-II Cómparacion I l estadística II--------...j 
/
Hot spots y Cold ~---L--+------+---------l1 índice de I 
spots de cambio I coincidencia difusa 
Resultato de la 
evaluación 
~/ Datos I 
1’955,577 habitantes, de los cuales el 84% se encontraba en localidades urbanas y el%16
en localidades rurales.
Michoacán de Ocampo: Localizado en la región central de México. Contaba para el año
2010 con una población de 4’584,471 habitantes, el 69% de ellos concentrados en zonas
urbanas y un 31% en zonas rurales [INEGI, 2010]. El Estado tiene un área total de 58,599
km2, el principal uso del suelo es el agŕıcola alcanzado a cubrir un 29% del área total. La
segunda cobertura con mayor extensión es la vegetación de selva baja caducifolia con un
28%, seguida de Bosque de pino encino con un 13%, el restante 30% se divide en otros
usos y coberturas [INEGI, 2013].
Jalisco: Es el estado con mayor población de los seleccionados, contando en el 2010 con
7’844.830 habitantes distribuidos un 87% en áreas urbanas y 13% en áreas rurales [INEGI,
2010]. El principal uso del suelo es el agŕıcola ocupando un 28% de los 78,588 km2 con los
que cuenta el estado. En segundo lugar se encuentra la vegetación de selva baja caducifolia
con un 19%, en el tercer lugar se encuentran las coberturas de bosque de encino y bosques
de encino pino, ambas con 14% del área total. El restante 25% se encuentra distribuido
en 33 coberturas y usos del suelo adicionales a los ya mencionados [INEGI, 2013].
3.2. Datos y Sofware
Los datos empleados para este estudio se resumen en el Cuadro 3-1. Los detalles de las
fuentes de datos se encuentran en el Anexo 1, el cual contiene las direcciones url que correspon-
den. Los datos en formato raster fueron remuestreados a una resolución de 250m. En el caso
de los datos continuos se empleó una interpolación bilineal. Para los datos categóricos se usó el
método del vecino más cercano.
El software empleado para gran parte de este estudio fue R [R Development Core Team,
2008], especialmente los paquetes Bfast [Verbesselt et al., 2012a,b] y Bfast spatial [Dutrieux,
2016]. Adicionalmente se emplearon Dinamica Ego 4.0, Postgis 2.1 y ArcGis 10.5.
24
Figura 3-2: Área de estudio
3.3. Imagenes MODIS - NDVI
El sensor MODIS 1, abordo de los satelites TERRA y AQUA, es un instrumento que cuenta
con 36 bandas espectrales, que van desde los 0.4 µm a los 14.4µm y con una resolución ra-
diométrica de 12 bits. Estas caracteŕısticas han permitido al equipo cient́ıfico detrás del proyec-
to generar diferentes productos con diferentes niveles de procesamiento, como los denominados
productos de valor agregado derivados de variables geof́ısicas mapeadas [Yengoh et al., 2015b].
Entre estos productos se encuentra el MOD13Q1 el cual contiene la capa de NDV I (́Indice de
Vegetación de Diferencia Normalizada, por sus siglas en inglés). El NDV I es posiblemente el
ı́ndice vegetación más usado, siendo éste la proporción de diferencia entre la banda del infra-
rrojo cercano -NIR (por sus siglas en ingles) y la banda roja-R, y la suma de estas dos bandas
(Ecuación 3-1) . [Chuvieco Salinero, 2010; Yengoh et al., 2015a]).
1
Sensor espectroradiómetro para imágenes de resolución moderada, MODIS por sus siglas en ingles.
25
Cuadro 3-1: Datos
No Dato de entrada Tipo Fuente Descripción
1 MOD13Q1-NDVI Raster NASA Composición de 16 d́ıas de NDVI -
250m
2 MOD13Q1-QA Raster NASA Capa de calidad del NDVI - 250m
3 MDE Raster INEGI Modelo Digital de Elevación - 90m
4 Cubierta arbórea Raster Hansen et al. [2013] Capa de cobertura arbórea global -
30m
5 Pérdida de bos-
que
Raster Hansen et al. [2013] Capa de eventos de pérdida de bos-
que (2000 - 2013) - 30m
6 Ganancia de bos-
que
Raster Hansen et al. [2013] Capa de eventos de ganancia de bos-
que (2000 - 2012) - 30m
7 CUS Vector INEGI Uso de suelo y vegetación serie V -
1:250,000
8 Localidades Vector INEGI Localidades urbanas y rurales
9 Red vial Vector INEGI Datos vectoriales de carreteras y
vialidades 1:50,000
10 Ŕıos Vector INEGI Red hidrográfica escala 1:50,000
11 Crecimiento Tabla [Ghilardi, 2008] Tasas de crecimiento asociadas a las
coberturas del mapa de CUS.
NDV I =
NIR�RED
NIR+RED
(3-1)
La capacidad de este ı́ndice se deriva del comportamiento radiométrico de la vegetación sana,
la cual refleja en la banda roja del espectro (0,6 a 0,7 µm) y en la banda del infrarrojo cercano
(0.7 a 1.1 µm).
El MOD13Q1 es la composición de los mejores valores observados por pixel en un periodo
de 16 d́ıas, con un tamaño de pixel de 250m. El algoritmo empleado considera que los pixe-
les contaminados por nubes y tomados cuando el sensor se encuentra fuera del nadir, tienen
menor calidad. Por tanto, un pixel libre de nubes, tomado sobre el nadir y sin contaminación
atmosférica residual es un pixel con buena calidad [NASA]. El algoritmo implementado solo
conserva los datos filtrados de buena calidad, sin embargo no todas las anomaĺıas son filtradas
por este algoritmo motivo por el cual existe una capa de calidad, que resume estas anomaĺıas
(QA-DSDs).
26
3.4. Preprocesamieto
3.4.1. Filtros de calidad
El preprocesamiento de las imágenes MODIS se realizó siguiendo metodoloǵıas empleadas
en estudios previos realizados en regiones como Sudamérica, Mongolia y China [Chuvieco, 2008;
Julien and Sobrino, 2010; Geng et al., 2014]. El primer paso del preprocesamiento consistió en
filtrar los datos de cada imagen, usando la información contenida en la capa QA-DSDs. Esta
capa contiene 16 bits, que describen la calidad de cada pixel,como se puede observar en el
Cuadro 3-2.
Los bits de 0 a 5 describen los parámetros generales de la calidad de cada ṕıxel, y en los
bits de 6 a 15 se describen de manera más especifica. Para este estudio se empleó un filtro que
considera la calidad general, utilizando una mezcla de nubes, máscara de tierra/agua y sombras.
En la Figura 3-3 se detalla el procedimiento utilizado.
Cuadro 3-2: Descripción de la capa de calidad de los ı́ndices de vegetación. Fuente:
http://www.ctahr.hawaii.edu/grem/mod13ug/sect0005.htm.
Bits Parámetro Valor Descripción
0-1 Calidad general 00 Producido con alta calidad
01 Producido con alta calidad revisar otros aspectos
10 Producido pero probablemente con nubes
11 Ṕıxel no producido por razones distintas a las nubes
2-5 Utilidad 0000 Perfecto
0001 Alta calidad
0010 Buena calidad
0011 Calidad aceptable
0100 Calidad aceptable mı́nima
0101 Calidad intermedia
0110 Inferior a calidad intermedia
0111 Calidad promedio
27
Cuadro 3-2: Descripción de la capa de calidad de los ı́ndices de vegetación. Fuente:
http://www.ctahr.hawaii.edu/grem/mod13ug/sect0005.htm.
Bits Parámetro Valor Descripción
1000 Inferior a la calidad promedio
1001 Calidad dudosa
1010 Calidad superior a datos excluidos
1011 Datos excluidos
1100 Baja calidad
1101 Sin corrección atmosférica
1110 Calidad demasiado baja para ser útil
1111 No es útil por otras razones, no se generó.
6-7 Cantidad de aerosoles 00 Clima
01 Bajo
10 Intermedio
11 Alto
8 Detección de nubes adyacentes 0 No
1 Si
9 Corrección atmosférica BRDF 0 No
1 Si
10 Mezcla de nubes 0 No
1 Si
11-13 Máscara tierra/agua 000 Océano superficial
001 Solo tierra
010 Costa de océanos y orilla de lagos
011 Aguas continentales superficiales
100 Aguas ef́ımeras
101 Aguas profundas continentales
110 Océano próximo a la costa
28
Cuadro 3-2: Descripción de la capa de calidad de los ı́ndices de vegetación. Fuente:
http://www.ctahr.hawaii.edu/grem/mod13ug/sect0005.htm.
Bits Parámetro Valor Descripción
111 Océano profundo
14 Posibilidad de hielo o nieve 0 No
1 Si
15 Posibilidad de sombras 0 No
1 Si
Figura 3-3: Algoritmo de filtrado empleando la capa de calidad de NDVI.
3.4.2. Interpolación y Filtro Savitzky Golay
Una vez filtrados los datos de cada ṕıxel, encontramos vaćıos en las series de tiempo, produc-
to de los valores que no alcanzan un grado de calidad suficiente. Por este motivo se completaron
los valores faltantes utilizando una interpolación ”spline”, esto soluciona el problema de vaćıos
29
en las series de tiempo, aunque continúa la contaminación de nubes, variabilidad atmosférica y
efectos bidireccionales [Chen et al., 2004]. Por esta razón se han formulado, aplicado y compa-
rado una serie de métodos orientados a la reducción del ruido, encontrando que el filtro para
señales Savitzky Golay presenta un gran desempeño [Chen et al., 2004; Julien and Sobrino,
2010; Wei et al., 2016].
El filtro Savitzky Golay se puede interpretar como un filtro de media móvil ponderada,
con una ponderación dada como un polinomio de cierto grado. Cuando se aplica a una señal,
se realiza un ajuste polinomial de mı́nimos cuadrados dentro de la ventana del filtro. Este
polinomio está diseñado para preservar momentos máximos dentro de los datos y para reducir
el sesgo introducido por el filtro. Este filtro se puede aplicar a cualquier dato consecutivo cuando
los puntos de los datos están en un intervalo fijo y uniforme a lo largo de la abscisa [Chen et al.,
2004]. Las series de tiempo de NDVI tienen estas propiedades.
La Ecuación 3-2 corresponde al filtro Savitzky Golay, donde Y es el valor original de NDVI,
Y ⇤ es el valor resultante de NDVI filtrado, Ci es el coeficiente para el i-ésimo valor NDVI
del filtro (ventana de suavizado), N es el numero de convoluciones y es igual a la ventana de
suavizado de tamaño (2m+1). El ı́ndice j es el ı́ndice de ejecución de la tabla de datos original.
Y ⇤j =
Pi=m
i=�mCiY
j+i
N
(3-2)
3.5. Bfast monitor
Esta investigación detectó los cambios sobre la cubierta vegetal empleando el algoritmo Bfast
monitor, el cual fue diseñado para el monitoreo de cambios en tiempo cuasi-real [Verbesselt et al.,
2012b,a]. Este algoritmo es una evolución del algoritmo Breaks For Additive Seasonaland Trend
- BFAST desarrollado por Verbesselt et al. [2010a,b], el cual tiene como propósito la detección
de perturbaciones en una serie temporal observada t = 1, ..., n. En contraste, Bfast monitor
busca si las nuevas observaciones t = n, n + 1, . . . mantienen el comportamiento esperado de
la muestra histórica t = 1, ...n. En otras palabras el método es capaz de identificar y modelar
30
Figura 3-4: Representación gráfica del periodo histórico y de monitoreo empleado para el Análi-
sis de Serie de Tiempo con Bfast.
automáticamente un historial estable, y luego detectar los cambios en los datos más recientes
Verbesselt et al. [2012b].
El intervalo de la serie de tiempo que se analiza para detectar los cambios abruptos se
denomina periodo de monitoreo. Tal como se mencionó anteriormente los cambios abruptos
se identifican en función del periodo histórico, el cual corresponde al intervalo que precede
al periodo de monitoreo. Sin embargo, en muchas ocasiones la totalidad del periodo histórico
puede no representar un comportamiento estable o normal. Por lo tanto se define un periodo
histórico representativo o estable como un intervalo del periodo histórico que esta libre de
cambios abruptos y se utiliza para detectar el cambio en el periodo de monitoreo.
Este estudio toma como periodo histórico el intervalo 2000 - 2007 y como periodo de moni-
toreo el intervalo 2008 - 2016 (Figura 3-4).
En resumen el BFAST monitor se puede explicar como un procedimiento de tres pasos
sucesivos:
Identificar un peŕıodo de histórico estable.
Modelar el peŕıodo histórico estable.
Analizar el periodo de monitoreo en búsqueda de cambios abruptos.
31
3.5.1. Modelo Tendencial y Estacional
Bfast monitor parte del modelo creado por [Verbesselt et al., 2012a,b], en el cual una serie
de tiempo es modelada mediante tres componentes, un componente estacional, un componente
tendencial y un término de error. 3-3
Yt = Tt+ St + "t (3-3)
Donde:
Tt = ↵1 + ↵2t (3-4)
La Ecuación 3-4 corresponde al componente tendencial, siendo ↵1 el intercepto y ↵2t la
pendiente de un modelo lineal a trozos.
St =
kX
j=1
� sin(
2⇡jt
f
+ �j) (3-5)
La Ecuación 3-5 corresponde a la función armónica que representa el componente estacional,
siendo �1, ..., �k la amplitud, �1, ..., �n la fase, f la frecuencia por unidad de tiempo (e.g f =23
observaciones anuales para una series de tiempo de 16 d́ıas) y k es el término de la función
armónica.
"t es un término de error que no puede ser observado en el tiempo t, es decir que no puede
ser explicado por ninguno de los componentes previamente mencionados.
Con la finalidad de ajustar el modelo a los datos del periodo histórico estable, la Ecuación
3-3 se reescribe como un modelo de regresión lineal [Verbesselt et al., 2012b].
yt = xTt � + "t,
xt = {1, t, sin(2⇡1tf ), cos(
2⇡1t
f ), ..., sin(
2⇡kt
f ), cos(
2⇡kt
f )}
T
� = {↵1,↵2, �1 cos(�), �1 sin(�), ..., �k cos(�), �k sin(�)}T
(3-6)
Donde � puede ser estimado usando mı́nimos cuadrados ordinarios.
32
3.5.2. Monitoreo del cambio
Empleando las ecuaciones vistas anteriormente es posible estimar un modelo de tendencia y
estacionalidad estable en un peŕıodo de tiempo observado. Ahora el punto central es determinar
si esta estabilidad se mantiene para nuevas observaciones [Verbesselt et al., 2012b]. Con esta
finalidad y dado que el modelo se puede formular como una regresión lineal (ver Ecuación 3-6),
los desarrolladores integraron técnicas de monitoreo de cambios estructurales empleando sumas
móviles acumuladas (MOSUM - por sus siglas en ingles.) sobre los residuos del periodo de
monitoreo t = n +1, .., n. Cuando el resultado de la Ecuación 3-7 se aleja de cero más allá de
un nivel de significancia estad́ıstica del 95%, Bfast monitor identifica éste como un punto de
ruptura.
MOt =
1
�̂
p
n
tX
s=t�h+1
(Ys � Ŷs) (3-7)
Donde �̂ es el estimado de varianza, Ys la observación actual, Ŷs la observación esperada,
n el número de observaciones en el periodo histórico y h el ancho de banda MOSUM definido
como fracción de n.
Independientemente de si detecta o no cambio Bfast monitor calcula la magnitud de cambio,
la cual se define como la mediana entre los valores esperados y los observados en el periodo de
monitoreo.
3.5.3. Peŕıodo de histórico estable
El éxito del enfoque de monitoreo descrito radica en que el periódico histórico t = 1, ..., n,
se encuentra libre de perturbaciones y puede ser empleado para modelar el comportamiento
normal esperado, empleando la Ecuación 3-3. En la practica esto ocurre con poca frecuencia,
por lo que en muchos casos no se usarán todas las observaciones disponibles y se empleará un
periodo histórico estable, inmediatamente anterior al inicio del periodo de monitoreo.
BFAST monitor permite identificar un peŕıodo histórico estable de manera manual o au-
tomática. La selección manual se puede realizar eligiendo un intervalo de tiempo antes del
periodo de monitoreo, de esta manera los datos desde el momento elegido hasta el comienzo del
33
periodo de monitoreo, se utilizan como un peŕıodo de historia estable. Este procedimiento se re-
comienda exclusivamente cuando se cuenta con conocimiento experto [Verbesselt et al., 2012a].
En los casos en los que no se cuenta con un conocimiento a priori sobre el comportamiento de
la variable objeto de estudio, o se cuenta con un gran numero de series de tiempo (como es
el caso de la imágenes de satélite), Bfast monitor tiene la capacidad de identificar un periodo
histórico estable automáticamente con un enfoque basado en datos.
Para definir de manera automáticamente el periodo de monitoreo se emplea una suma
acumulada de residuos en orden inverso (ROC-por sus siglas en inglés), que evalúa el error en
la predicción de la ecuación 3-3 al retroceder en el tiempo t = n, n � 1, n � 2, . . . hasta que se
presenta una ruptura.
La Figura 3-5 es el resultado del algoritmo Bfast monitor para una serie de tiempo extráıda
del NDVI del productor MOD13Q1. En esta gráfica se puede observar los datos del periodo
histórico (linea negra continua), el modelo ajustado en función del periodo histórico estable
(linea azul), el inicio del periodo de monitoreo (linea negra discontinua), los nuevos datos para
el periodo de monitoreo (linea roja continua) y un punto quiebre detectado en el mes de abril
de año 2011.
3.5.4. Bfast Spatial
Bfast Spatial es una herramienta desarrollada en el año 2016 [Dutrieux, 2016], la cual aplica
el enfoque de Bfast monitor, en un contexto de imágenes de satélite, donde en potencia cada
ṕıxel corresponde a una serie de tiempo.
El resultado de BFAST Spatial es una pila de tres capas: 1) Los puntos de ruptura detectados
para cada ṕıxel, con información sobre la temporalidad de ocurrencia en formato de año decimal;
2) La magnitud de cambio para cada ṕıxel; y 3) Un indicador de error - posibles valores de 1 si
se ha encontrado un error en un ṕıxel particular o NA donde el algoritmo tuvo éxito.
Para detectar los cambios sobre la cubierta vegetal se consideran sólo la magnitud de los
ṕıxeles donde ocurrió una ruptura, las magnitudes negativas están asociadas con procesos defo-
34
Figura 3-5: Ejemplo de la salida gráfica del algoritmo Bfast monitor aplicado en una zona de
selva baja subcaducifolia en Michoacán.
restación y degradación, mientras que las magnitudes positivas están asociadas a revegetación
o reforestación. Sin embargo la presencia o ausencia de puntos de ruptura para definir los cam-
bios no es suficiente, ya que los puntos de ruptura asociados a magnitudes positivas o negativas,
pero los próximos a cero no se traducen directamente en cambios [Dutrieux et al., 2016, 2015].
Para evitar la detección de falsos positivos es necesario establecer un umbral. Este estudio
empleará un umbral de 0.05 (unidades de NDVI) de acuerdo a lo sugerido por Dutrieux [2016].
Los resultado obtenidos en 3.4.2, fueron procesados con Bfast spatial, para encontrar los
cambios sobre la cubierta vegetal.
3.6. Modelización de escenarios de ahorro de leña - MoFuSS
MoFuSS es un modelo dinámico que tiene como objetivo simular el efecto espacio temporal
de la recolección de leña sobre la vegetación local, y que representa el ahorro de biomasa no
renovable debido a una reducción en el consumo por alteración en los patrones de cosecha y
35
recrecimiento en el tiempo [Ghilardi et al., 2016].
MoFuSS está compuesto por cuatro módulos (ver Figura 3-6): 1) un componente de fricción
que crea mapas de impedancia; 2) un algoritmo modificado de distancia inversa ponderada
(IDW - por sus siglas en inglés), que crea mapas de presión sobre la propensión a los eventos
de cosecha de leña; 3) un componente de oferta/demanda que proyecta la cantidad esperada
de leña que se va a cosechar en cada peŕıodo de tiempo para cada ṕıxel, y la respuesta de la
vegetación a esa perturbación; y 4) un módulo de pérdida y ganancia de bosque que proyecta
el desmonte o la ganancia en eventos esperados.
Figura 3-6: Algoritmo MoFuSS simplificado, En el cual se observan las entradas y procesos
principales dentro de cada uno de los módulos, y los resultados finales. Fuente: Ghilardi et al.
[2016]
3.6.1. Mapas de fricción
La fricción o impedancia fue abordado inicialmente dentro de las ciencias de la información
geográfica para realizar modelos de redes, en donde ésta es una medida de la resistencia o
coste, requerida para recorrer una ruta en una red. Esta resistencia puede ser una medida de
36
la distancia de viaje, el tiempo, la velocidad de viaje multiplicado por la distancia, etc., por
mencionar algunas [Djokic and Maidment, 1993]. Los valores de impedancia más altos indican
más resistencia al movimiento y un valor de cero indica que no hay resistencia.
Para MoFuSS el mapa de fricción caracteriza en cada ṕıxel el tiempo que puede necesitar
un recolector de leña para alcanzar un punto de cosecha a pie o en automóvil (dependiendo de
las rutas utilizadas en una localidad determinada). Los datos requeridos para construir estos
mapas de fricción son las velocidades de desplazamiento de los recolectores de leña para caminar
y conducir, afectados por caracteŕısticas topográficas tales como la pendiente, las condiciones
de las carreteras, los ŕıos y los cuerpos de agua, o los tipos de cubierta terrestre.
3.6.2. Componente IDW
El componente IDW crea un mapa de presión, representando la probabilidad o propensión
de cada ṕıxel a ser cosechado para leña, para dos tipos de colectores de leña: 1) personas que
viajan a pie y recogen leña para uso doméstico, y 2) vendedores de leña comercial que usan
veh́ıculos, lo que les permite acceder a áreas distantes y transportar grandes volúmenes de leña
(ver Ecuación 3-8).
P(t)k =
⇣ nX
i=1
Cik
dnik
⌘
(t)
(3-8)
Donde P(t)k es un ı́ndice de la presión por la colecta de leña k para cualquier momento
en el tiempo t, C es el consumo de leña residencial en tDM por localidad; y n es número
positivo real que controla la función de decadencia de la interpolación, en otras palabras define
el peso de los puntos vecinos. La gran ventaja de este ı́ndice es que cada ṕıxel dentro del
área de estudio esta influenciado por los centros de consumo, es decir que a cada ṕıxel que
aloja cualquier cantidad de leña cosechable se le asigna una probabilidad de cosecha, que es
directamente proporcional a la demanda de leña de todas las localidades dentro del área de
estudio, e inversamente proporcional a la proximidad basada en la fricción.
37Fuente Descripción
Bosques Leña que se obtiene mediante la poda o
tala de árboles vivos o la recolección de
ramas cáıdas y madera muerta.
Árboles fuera de los bosques
(TOF)
Incluye árboles en tierras de cultivo, com-
puestos domésticos y bienes comunes al
borde de las carreteras, a los que se acce-
de mediante la poda de árboles vivos y/o
la recolección de ramas muertas o descen-
dientes.
Actividades de desmonte Incluye la tala de bosques o matorrales
para nuevos cultivos o pastos y constitu-
ye una importante fuente de suministro de
madera.
Cuadro 3-3: Descripcion de las fuentes de oferta de leña.
3.6.3. Demanda y Oferta de leña
La demanda de leña es calculada por MoFuSS empleando la Ecuación 3-9.
Ct = (
nX
i=1
nX
i=j
hhijujfci)t + (
nX
i=1
nX
i=j
hhijujfbi)t (3-9)
Donde Ct es el total de leña residencial en toneladas de masa seca (tDM - por sus siglas en
inglés) para cualquier momento en el tiempo t, hhij es el numero de hogares que usan leña por
comunidad i, utilizando un dispositivo de cocción j, uj es el consumo de los hogares en tDM,
fci y fbi el promedio de leña cosechada y comprada, respectivamente.
En cuanto a la oferta de leña se consideran principalmente de tres fuentes, Bosques, Árboles
fuera de los bosques (TOF - por sus siglas en inglés) y Actividades de desmonte (ver Cuadro
3-3) .
La localización de estas fuentes se obtienen de un mapa de cobertura y uso del suelo, el cual
es uno de los insumos requeridos por MoFuSS.
38
3.6.4. Pérdida y ganancia de bosque.
MoFuSS supone que el crecimiento de la biomasa leñosa se da en función de tres factores:
1) El almacén de biomasa en el estado anterior, 2) la tasa de crecimiento máxima (rmax) y 3) la
densidad máxima de biomasa o capacidad de carga (K) [Bailis et al., 2015]. rmax y K dependen
directamente de factores biof́ısicos como el tipo de suelo, la precipitación, la insolación y la
altitud, entre otros. Estos parámetros raramente están disponibles, por lo que MoFuSS recibe
los parámetros rmax y K como una tabla asociada a las clases del mapa de cobertura y uso del
suelo. En la Ecuación 3-10 se puede observar la relación entre los parámetros mencionados.
AGB(t+1)i = AGB(t)i +AGB(t)irmax
⇣
1�
AGB(t)i
Ki
⌘
(3-10)
Donde AGB(t)i
2 y AGB(t+1)i corresponden a la biomasa leñosa sobre el suelo disponible
para leña en cada clase i del mapa de cobertura y uso del suelo, en el tiempo t y t + 1,
respectivamente.
Para AGBt=0 se recomienda emplear mapas que contengan una distribución espacial conti-
nua de biomasa tales como las producidas por Avitabile et al. [2016]; Baccini et al. [2015], ya
que estos ofrecen mayor precisión y se encuentran calculados por ṕıxel.
Adicionalmente un ṕıxel puede volverse incosechable si la AGB cae por debajo del umbral
determinado por el modelador, hasta que su tasa de recrecimiento lo lleve a superar este um-
bral. Para definir dicho umbral en la actualidad se carece de estudios detallados, por lo que
usualmente es definido mediante conocimiento experto. Ghilardi et al. [2016] Se recomienda un
umbral de 5.0 tDMha para colectores tipo 2 y 0.1
tDM
ha para colectores tipo 1.
MoFuSS emplea simulaciones de Monte Carlo (MC), para introducir la incertidumbre aso-
ciada a los parámetros de crecimiento de la AGB. Con cada iteración de MC los parámetros de
crecimiento de AGB vaŕıan de manera aleatoria. Las simulaciones de MC se emplean también
sobre la distribución espacial de ṕıxeles visitados por los menos una vez, y un factor de corte
que regula el grado de aleatoriedad del mecanismo de siembra, el cual es empleado por un
2
AGB - biomasa aérea por sus siglas en inglés
39
algoritmo genético para simular la cosecha.
En esencia el modelo genera tres resultados: 1) la AGB restante (crecimiento menos cosecha,
para t = n), 2) NRB calculado por ṕıxel donde se han producido disminuciones en la AGB y
3) fNRB fracción del consumo total de leña no renovable. Como se mencionó en el capitulo 1
este estudio tomará el NRB calculado para la validación del modelo MoFuSS.
3.7. Análisis de puntos calientes y puntos fŕıos
El análisis de puntos calientes y puntos fŕıos se realiza mediante el cálculo de la medida
de asociación espacial Getis-Ord (ecuación 3-11) Gi⇤ [Getis and Ord, 1992; Ord and Getis,
1995; Chun and Gri�th, 2013], cuya idea fundamental es la de encontrar clusters espacialmente
expĺıcitos, buscando valores altos que se encuentren rodeados de valores altos (puntos calientes)
y valores bajos rodeados de valores bajos, en el caso de los puntos fŕıos.
Gi⇤ =
Pn
j=1wijxj � X̄
Pn
j=1wij
S
s⇥
nw2ij �
�Pn
j=1wij
�2⇤
n� 1
(3-11)
Donde S es la desviación estándar de la variable x y w es el peso de que tiene la localiza-
ción de i sobre la variable en la localización j [Chun and Gri�th, 2013]. Los valores positivos
estad́ısticamente significativos, indican que cuanto mayor es el valor Gi⇤ más fuerte será el clus-
ter de valores altos (puntos calientes). Para los valores negativos estad́ısticamente significativos,
cuanto menor sea el valor, será mas fuerte el cluster de los valores bajos (puntos fŕıos) [Chun
and Gri�th, 2013; Reddy et al., 2016].
Aplicar este análisis a las capas de pérdida de NDVI tiene como finalidad reclasificar esta
información de forma categórica, y de esta manera facilitar la evaluación del desempeño del
modelo MoFuSS.
40
3.8. Validación
La validacíıon de MoFuSS se realizo siguiendo tres estrategias la primera de ellas consiste
en aplicar análisis de correlaciones globales y locales, la segunda estrategia es una comparación
estad́ıstica empleando pruebas de hipótesis; finalmente la tercera estrategia es una comparación
espacial empleando un ı́ndice de coincidencia difusa.
3.8.1. Análisis de correlación
El análisis de correlación es el primer análisis aplicado con la finalidad de comparar los
resultados de pérdida de NDVI obtenidos del algoritmo Bfast Spatial y MoFuSS (3.5.4 y 3.6.4
respectivamente). El resultado esperado este análisis es conocer el grado o la intensidad de la
asociación [Steel et al., 1985] entre las capas mencionadas anteriormente [Steel et al., 1985]. Las
pruebas empleadas fueron las Pearson y Spearman, de manera global y usando ventanas móviles
de diferentes tamaños. Estos análisis de correlación se emplearon como análisis exploratorios ya
que las evaluaciones de modelos espaciales ṕıxel a ṕıxel son poco recomendadas [Hagen, 2003;
van Vliet et al., 2016; Koch et al., 2015; Homolová et al., 2014; Mas et al., 2012].
3.8.2. Contraste de hipótesis
Las pruebas de contraste de hipótesis son pruebas comunes en el análisis estad́ıstico, estas
tienen como propósito comparar muestras aleatorias de dos poblaciones con finalidad de deter-
minar si las poblaciones son estad́ısticamente iguales. En los casos en los cuales las muestras
de ajustan a una distribución normal se emplean pruebas como la t-student para comparar las
medias, la hipótesis nula para esta prueba indica que H0 : µ1 = µ2 y como hipótesis alternativa
se tiene que HA : µ1 6= µ2.
Cuando la muestra o la población a analizar no presenta una distribución normal se recurre
entonces a pruebas no parametricas tales como el estad́ıstico kruskal wallis o la U de Mann-
Whitney. Estas pruebas en esencia buscan lo mismo que las pruebas paramétricas, pero modifi-
can la hipótesis realizando la comparación sobre la mediana y no sobre la media. Adicionalmente
41
las pruebas no paramétricas son menos potentes que las pruebas paramétricas [Burt, 2009].
3.8.3. Índice de coincidencia difusa
Los análisis de similitud entre capas raster categóricas, tradicionalmente se generaban rea-
lizando análisis de coincidencia ṕıxel por ṕıxel. Este método puede registrar una discrepancia
entre ṕıxeles, incluso ante desplazamientos menores entre los ṕıxeles y donde los mapas ca-
tegóricos representan en esencia un mismo patrón espacial[Yan and Li, 2015]. Para resolver
este problema Hagen [2003], desarrolló una métrica de coincidencia difusa (fuzzy), la cual ha
sido empleada ampliamente en la evaluación de modelos espacialmente expĺıcitos, como los mo-
delos de cambio de cobertura y uso del suelo, hidrológicos y de servicios ecosistémicos [van Vliet
et al., 2016; Koch et al., 2015; Homolová et al., 2014; Mas et al., 2012].
El objetivo de esta métrica es obtener un análisis espacial y gradual de la similitud entre dos
raster categóricos dentro de un contexto de vecindad. Esta prueba se fundamenta en el concepto
de difuminación (fuzzines) espacial, donde la representación de un ṕıxel es influenciada por el
propio ṕıxel y por los ṕıxeles de su vecindad [Mas et al., 2012; Hagen, 2003]. Una prueba de
coincidencia difusa bidireccional está implementada en el software DINAMICA Ego, aplicando
teoŕıa de conjuntos difusos y una función de decaimiento constante dentro de un tamaño de
ventana variable. Si el mismo número de celdas de una categoŕıa determinada se encuentra
dentro de la ventana, el ajuste será 1 independientemente de su ubicación. Esto representa una
manera conveniente de evaluar la aptitud del modelo a través de la disminución de la resolución
espacial. Los modelos que no coinciden bien en alta resolución pueden tener un ajuste apropiado
a una resolución más baja.
42
Caṕıtulo 4
Resultados y discusión
4.1. Estimación de pérdida de la cobertura vegetal
La Figuras 4-1, 4-2 y 4-3 muestran los resultados del algoritmo Bfast monitor (Ecuación 3-3)
aplicada a las series de tiempo para los estados de Yucatán, Michoacán y Jalisco. En estas figuras
se pueden observar las magnitudes de NDVI negativas menores a -0.05, que coinciden con ṕıxeles
en los cuales ocurrió una ruptura en las series de tiempo, representando una pérdida sobre la
cobertura para el peŕıodo 2008 - 2016. Adicionalmente con la finalidad de evitar sobrestimar
los cambios, se aplicó una máscara que excluye las coberturas no susceptibles a la cosecha de
leña (e.g. agricultura, plantaciones forestales, áreas urbanas).
4.2. Impactos proyectados por extracción de leña a partir de
MoFuSS
Las simulaciones en MoFuSS se realizaron para el periodo 2008 - 2016. En las Figuras 4-4,
4-5 y 4-6 se observan los valores de NRB expresados en toneladas por ṕıxel. Los resultados
adicionales del modelo se pueden consultar el en Anexo 2, el cual contiene los informes generados
por MoFuSS para cada una de las tres simulaciones (estados de Yucatán, Michoacán y Jalisco).
43
Figura 4-1: Pérdida sobre la cubierta vegetal en unidades de NDVI para el estado de Yucatán
(2008 - 2016)
Figura 4-2: Pérdida sobre la cubierta vegetal en unidades de NDVI para el estado de Michoacán
(2008 - 2016)
44
Figura 4-3: Pérdida sobre la cubierta vegetal en unidades de NDVI para el estado de Jalisco
(2008 - 2016)
Figura 4-4: Estimacion de la extracción de biomasa leñosa a tasas que exceden la regeneración
natural (NRB) dentro del peŕıodo 2008 - 2016, para el estado de Yucatán.
45
Figura 4-5: Estimación de la extracción de biomasa leñosa a tasas que exceden la regeneración
natural (NRB) dentro del peŕıodo 2008 - 2016, para el estado de Michoacán.
Figura 4-6: Estimación de la extracción de biomasa leñosa a tasas que exceden la regeneración
natural (NRB) dentro del peŕıodo 2008 - 2016, para el estado de Jalisco.
46
4.3. Agrupaciones de áreas con pérdida de NDVI (Puntos ca-
lientes - puntos fŕıos)
Previamente al cálculo de los clusters (puntos calientes y fŕıos) de pérdida se aplicó un cálculo
de valor absoluto con la finalidad de facilitar posteriormente la interpretación de los resultados.
La Figuras 4-7, 4-8 y 4-9 muestran los resultados del ı́ndice Getis-Ord con un p < 0.05.
Figura 4-7: Puntos fŕıos y puntos calientes de pérdida de NDVI, calculados mediante el ı́ndice
Getis-Ord para el estado de Yucatán
47
Figura 4-8: Puntos fŕıos y puntos calientes de pérdida de NDVI, calculados mediante el ı́ndice
Getis-Ord para el estado de Michoacán
Figura 4-9: Puntos fŕıos y puntos calientes de pérdida de NDVI, calculados mediante el ı́ndice
Getis-Ord para el estado de Jalisco
48
4.4. Resultados de las estrategias de validación de MoFuSS
4.4.1. Análisis de correlación entre Biomasa No Renovable (NRB) y pérdida
del NDVI
El primer paso de la validación fue la comparación entre los valores de NRB obtenidos de
MoFuSS y los valores de pérdida de NDVI. Este paso consistió en buscar correlaciones globales
y locales entre estos dos raster. En el caso de las correlaciones globales los coeficientes de
correlación son cercanos a 0 (ver Cuadro 4-1). En cuanto a las correlaciones locales éstas se
pueden observar en la Figura 4-10, y se emplearon en ellas cuatro diferentes tamaños para las
ventanas móviles.
Estado Coeficiente
Yucatán 0.06
Michoacán 0.02
Jalisco 0.13
Cuadro 4-1: Correlaciones globales entre los valores de NRB y la perdida de NDVI.
4.4.2. Comparación estad́ıstica
Se realizó un análisis de autocorrelación espacial de NRB, con la finalidad de tomar mues-
tras no autocorrrelacionadas espacialmente, para posteriormente comparar los valores de NRB
al interior de los clusters de pérdida con los valores fuera de éstos. Los resultados de los semi-
variogramas ajustados a un modelo esférico se pueden observar en el Cuadro 4-2, los resultados
obtenidos indican una alta autocorrrelación espacial por lo cual no es posible obtener una
muestra no autocorrelacionada. Ante esto se decidió emplear una muestra de 15,000 ṕıxeles de
NRB.
Estado Nugget Modelo Sill Rango (m)
Yucatán 3.749622 Sph 2.802967 81,739.73
Michoacán 1.494724 Sph 3.513463 110,795.00
Jalisco 0.2443015 Shp 0.3095379 109,238.40
Cuadro 4-2: Autocorrelación espacial.
49
F
ig
u
ra
4-
10
:
C
or
re
la
ci
on
es
lo
ca
le
s
en
tr
e
lo
s
va
lo
re
s
d
e
N
R
B
y
la
p
ér
d
id
a
d
e
N
D
V
I.
V
en
ta
n
as
m
óv
il
es
d
e
a)
3x
3
ṕ
ıx
el
es
,
b
)
7x
7,
c)
9x
9
y
d
)
11
x
11
50
e 
'" 
'ro 
"O U 
::> ro 
o .~ -I:: 
.3 .~ 
~ 
~ 
§ 
~ <i " ci -: 
~ 
~ 
ª 
z z z z z 
f<l 
:¡, ~ :¡, e 
~ ~ 
apn¡I¡Bl 
"' ~ 
~ 8 
.~ ~ 
~ 
o ro 
...J 2-
o 
~ ~ <i 
o 
~ -: 
~ 
~ 
;: 
~ 
z 
~ 
z z 
1::: f<l 
apn¡I¡Bl 
'" e "O '2 ::> 
.~ ro u 
.3 ::J e 
;: 
iil 
~ 
::i 
~ 
~ 
¡¡i 
~ 
z z z z z 
:¡, 
N 
:¡, :¡, 
N ¡;j ~ 
apn¡!¡el 
La prueba de normalidad Kolmogorov–Smirnov arrojó que la distribución de NRB para
ninguno de los estados se ajusta a una distribución normal (ver Cuadro 4-3). Los grupos a
comparar son los valores de NRBpf (valores de NRB en el cluster de puntos fŕıos de pérdida),
NRBpc (valores de NRB en el cluster de puntos calientes) y los valores de NRB que se en-
cuentran fuera de estos clusters. Una comparación descriptiva de estos grupos compuesto por
5,000 ṕıxeles, se puede observar en los diagramas de cajas de la Figura 4-11. Adicionalmente
en el Cuadro 4-4 se encuentran las estad́ısticas descriptivas para cada uno de ellos.
La comparación entre los grupos se realizó inicialmente empleando la prueba Kruskal-Wallis
la cual tiene como hipótesis nula H0 : ⌘1 = ⌘2. Los resultados de esta prueba que se presentan
en el cuadro 4-5 indican que existe diferencia en al menos 2 de los 3 grupos comparados. Una
prueba U de Mann-Whitney fue empleada para comprobar la diferencia entre los clusters y los
valores de NRB fuera de estos (Cuadro 4-6). Se encontró en los tres estados analizados, que
existen diferencias estad́ısticamente significativas entre NRB fuera de los clusters y los valores
de NRB al interior de los clusters de puntos fŕıos (NRBpf).
Estado Estad́ıstico Valor de p
Yucatán 0.14688 < 2.2e� 16
Michoacán 0.14754 < 2.2e� 16
Jalisco 0.18922 < 2.2e� 16
Cuadro 4-3: Prueba de normalidad a los valores de NRB.
Estado Grupo Media Desviación estándar Mediana
Yucatán NBR 1.91 1.82 1.35
NRBpc 1.87 1.82 1.27
NRBpf 3.38 2.79 2.75Michoacán NBR 5.61 3.42 5.57
NRBpc 5.65 3.44 5.63
NRBpf 3.42 2.53 2.91
Jalisco NBR 1.06 0.89 0.78
NRBpc 1.05 0.9 0.74
NRBpf 0.75 0.74 0.44
Cuadro 4-4: Estad́ısticas descriptivas para cada uno de los grupos establecidos. NRBpf ) valores
de NRB al interior del cluster de punto fŕıos de perdida, NRBpc ) valores de NRB al interior
cluster de puntos calientes, NRB) valores de NRB fuera de los clusters.
51
(a) Yucatán (b) Michoacán
(c) Jalisco
Figura 4-11: Diagramas de cajas para cada uno de los grupos establecidos. NRBpf ) valores
de NRB al interior del cluster de punto fŕıos de perdida, NRBpc ) valores de NRB al interior
cluster de puntos calientes, NRB) valores de NRB fuera de los clusters.
52
Estado Estad́ıstico Grados de libertad Valor de p
Yucatán 1184 2 < 2.2e� 16
Michoacán 1397.5 2 < 2.2e� 16
Jalisco 494.77 2 < 2.2e� 16
Cuadro 4-5: Comparación entre los tres grupos definidos empleando la prueba de Kruskal-Wallis.
Estado Grupos Estad́ıstico Valor de p
Yucatán NRBpf - NRB 16741000 < 2.2e� 16
NRBpc - NRB 12331000 0.2408
Michoacán NRBpf - NRB 7877200 < 2.2e� 16
NRBpc - NRB 12637000, 0.3423
Jalisco NRBpf - NRB 9653200 1.821e� 08
NRBpc - NRB 12360000 0.3313
Cuadro 4-6: Comparación entre los valores de NRB al interior de los clusters y los valores fuera
de los cluesters, empleando la prueba U de Mann-Whitney.
4.4.3. Índice de coincidencia difusa entre áreas de NRB y Agrupaciones de
áreas con pérdida de NDVI
Mediante el ı́ndice de coincidencia difusa, se comparó la coincidencia espacial entre los
clusters de pérdida (PEDRDIDApf y PERDIDApc) contra las áreas de NRB. Adicionalmente
se generaron capas raster aleatoriamente que siguen la misma distribución estad́ıstica de NRB,
incluyendo máximos y mı́nimos. Estos datos aleatorios fueron comparados con los clusters para
evaluar si la coincidencia de espacial de las zonas de cosecha presentan una mejor coincidencia
que los datos aleatorios (Ver Figura 4-15).
Para evitar falsas coincidencias se emplearon diferentes umbrales para los valores de NRB
y se eliminaron agrupaciones de ṕıxeles inferiores 0.5 km2. Este procedimiento se aplicó sobre
las capas resultantes de MoFuSS (NRB) y sobre los valores generados aleatoriamente (NRBa);
encontrando el mejor desempeño para NRB > 2Tn en el caso de Yucatán y Michoacán. Estos
resultados se puede observar en las Figuras 4-16, 4-17. En el caso del estado de Jalisco no se
encontró un umbral en el cual el desempeño del modelo sea mejor que el desempeño de los datos
aleatorios como se puede observar en la Figura 4-18.
53
Figura 4-12: Comparación NRB generado por MofuSS (panel izquierdo) y NRB aleatorio (panel
derecho) para el estado de Yucatán
Figura 4-13: Comparación NRB generado por MofuSS (panel izquierdo) y NRB aleatorio (panel
derecho) para el estado de Michoacán
54
Figura 4-14: Comparación NRB generado por MofuSS (panel izquierdo) y NRB aleatorio (panel
derecho) para el estado de Jalisco
Figura 4-15: Esquema de la comparación espacial entre las áreas de NRB y NRB aleatorio con
las áreas definidas como clusters de perdida.
55
Figura 4-16: Relación entre el ı́ndice coincidencia difusa (eje Y) y las comparaciones entre el
NRB de MoFuSS y los clústers de pérdida de NDVI, y entre el NRB aleatorio (NRBa) y los
clústers de pérdida de NDVI; en función de la distancia (eje X), para el estado de Yucatán.
56
Figura 4-17: Relación entre el ı́ndice coincidencia difusa (eje Y) y las comparaciones entre el
NRB de MoFuSS y los clústers de pérdida de NDVI, y entre el NRB aleatorio (NRBa) y los
clústers de pérdida de NDVI; en función de la distancia (eje X), para el estado de Michoácan.
57
Figura 4-18: ÍRelación entre el ı́ndice coincidencia difusa (eje Y) y las comparaciones entre el
NRB de MoFuSS y los clústers de pérdida de NDVI, y entre el NRB aleatorio (NRBa) y los
clústers de pérdida de NDVI; en función de la distancia (eje X), para el estado de Jalisco.
58
4.5. Discusión
Los insumos usados para la evaluación de un modelo se asumen como los más cercanos a la
realidad posible. En este caso el análisis de series de tiempo permitió identificar zonas en las
cuales las variaciones en el NDVI indican perturbaciones forestales, que no se deben únicamente
a la cosecha de leña e involucra otros factores tales como la ganaderia, agricultura, la extracción
de madera con fines no energéticos, entre otros. Adicionalmente el cálculo de la disminución de
NDVI se realizó de manera genérica, es decir no se realizaron ajustes a la función armónica para
cada tipo vegetación. Estas consideraciones son importantes para discutir la validez del modelo,
ya que el umbral de cambio del NDVI podŕıa corresponder a un valor diferente de acuerdo al
tipo de vegetación.
4.5.1. Validación de MoFuSS usando series de tiempo de NDVI
Como se mencionó en el Capitulo 2 la validación del modelo se realizó siguiendo tres estrate-
gias. La primera de ellas fue el análisis de correlación entre los valores de NRB y los de pérdida
de NDVI. Los resultados eran hasta cierto punto esperados, ya que este tipo de evaluación ṕıxel
a ṕıxel ha demostrado ser poco eficaz en la validación de modelos espacialmente expĺıcitos [van
Vliet et al., 2016; Koch et al., 2015; Homolová et al., 2014; Mas et al., 2012; Hagen, 2003]. Este
análisis se realizó de forma exploratoria.
La segunda estrategia fue la comparación estad́ıstica de las medianas de NRB, mediante
pruebas de hipótesis. En éstas pruebas se demostró que existe una diferencia entre la respuesta
de la vegetación a la cosecha mediana de leña al interior de los clusters de puntos fŕıos, y
los valores de NRB fuera de éstos. Es decir, existe una asociación entre los valores de bajas
perturbaciones sobre la vegetación y los valores de cosecha predichos por el modelo.
Para los estados de Yucatacán y Michoacán las medianas indican una diferencia de 1.4tn y
2.66tn, entre los valores de NRB al interior de los clusters de puntos fŕıos y los valores fuera de
éstos. En el caso de Jalisco esta diferencia es 0.34tn, por lo que a pesar de ser estad́ısticamente
diferentes, en términos prácticos esta diferencia no es representativa.
59
Por otra parte los resultados de la comparación entre los valores de NRB al interior de los
clusters de puntos calientes y los valores de NRB fuera de éstos, demostraron que no existe
diferencia entre las medianas de estos grupos. Es decir que no existe una asociación entre los
lugares de mayor perturbación y las predicciones realizadas por el modelo.
La tercera y última estrategia fue la comparación espacial usando el ı́ndice de coinciden-
cia difusa. El resultado de esta comparación para el estado de Yucatán indica que el modelo
está generando predicciones con una alta coincidencia con las zonas de baja perturbación ob-
servadas en los ı́ndices de vegetación. Para el estado de Michoacán los resultados indican una
baja coincidencia entre los clusters de puntos fŕıos y las áreas de cosecha, sin embargo esta
coincidencia es mayor que la generada por el modelo aleatorio.
Esto concuerda con el conocimiento previo que se tiene sobre el impacto de la cosecha de
leña con fines residenciales en América Latina, es decir que este tipo de cosecha no es un factor
determinante de deforestación [Masera et al., 2015]. En cambio el rendimiento del modelo en
el estado de Jalisco tiene otro comportamiento ya que la coincidencia de NRBa, es mayor que
las estimaciones realizadas por el modelo. Esto se puede deber a lo siguiente:
Estudios anteriores han llegado a la conclusión de que en México la cosecha de leña tiene
una alta heterogeneidad espacial [Ghilardi et al., 2007], en otras palabras el patrón de
cosecha de leña puede variar ampliamente en función de la localización del fenómeno. De
acuerdo con esto puede ser que los patrones de cosecha en Jalisco generen impactos sobre
la cubierta vegetal a