jueves 12 de enero de 2012

Índices de Vegetación

Índices de Vegetación con Imágenes de Satélite 
Vegetation Index from Satellite Imagery

El análisis de la vegetación y la detección de los cambios en los patrones de vegetación son claves para la evaluación y el monitoreo de recursos naturales. Entonces no resulta sorpresivo que la detección y la evaluación cuantitativa de la vegetación verde es una de las aplicaciones principales de la percepción remota para el manejo de recursos ambientales y la toma de decisiones (Eastman 2003).

El empleo de los cocientes para discriminar masas de vegetales se deriva del peculiar comportamiento radiométrico de la vegetación. La signatura espectral característica de la vegetación sana muestra un claro comportamiento entre las bandas roja (0.6 a 0.7 µm) y el infrarrojo cercano (0.7 a 1.1 µm). Se produce  un notable contraste espectral entre la banda R del espectro y la del IRC, lo que permite separar la vegetación sana de otras cubiertas (Chuvieco 2008).

Figura1: Contraste espectral de la vegetación sana entre las bandas R e IRC del espectro

Cuando la vegetación sufre algún tipo de estrés, su reflectividad será inferior en el IRC, aumentando paralelamente en el rojo (al tener menor absorción clorofílica), con lo que el contraste en ambas capas será mucho menor. Los índices mas empleados son el cociente simple entre las bandas (Ci), y el denominado índice de vegetación de diferencia normalizada (NDVI) propuesto originalmente por Rouse et al (1974).

Índices de vegetación

Cociente simple (Ci)
El cociente simple (Ci) representa la relación entre las reflectividades del infrarrojo cercano y del rojo, los cuales representan las reflectividades de la banda 4 y 3 respectivamente, para el TM y ETM+ de las imágenes Landsat.
Índice de vegetación de diferencia normalizada (NDVI)
El Índice de Diferencia de Vegetación Normalizado, también conocido como Normalized Difference Vegetation Index (NDVI) (Rouse et al., 1974) por sus siglas en inglés. Es un índice usado para medir la diferencia normalizada entre las reflectancias del rojo y del infrarrojo cercano, proporcionando una medida sobre la cantidad, calidad y desarrollo de la cobertura vegetal y vigorosidad en áreas extensas.

Guyot y Gu, citados por Chuvieco (2008) usando un modelo teórico concluyen que los valores del NDVI para imágenes Landsat y SPOT calculados a partir de los ND subestiman entre 0.05 y 0.20 el valor calculado con reflectividades, siendo este error mayor con valores de NDVI inferiores a 0.5 y para las imágenes SPOT. En consecuencia, proponen una corrección que haga equivalente el cálculo con ND o reflectividades para las imágenes Landsat:

Índice de vegetación transformado (TVI)
Dering y colaboradores citados por Eastman (2003), modifica el NDVI agregando una constante de 0.5 a todos sus valores y calculando la raíz cuadrada de los resultados. La constante 0.5 se introduce para evitar operar con valores negativos del NDVI. El cálculo de la raíz cuadrada se emplea para corregir los valores del NDVI que se aproximan a una distribución Poison e introducir una distribución normal. Con estos dos elementos, el TVI toma la siguiente forma:

Índice de vegetación transformado corregido (CTVI)
Propuesto por Perry y Lautenschlager citados por Eastman (2003), apunta a corregir el TVI. Resulta obvio que agregar una constante de 0.5 a todos los valores del NDVI no siempre elimina los valores negativos porque los valores del NDVI pueden tener el rango -1 a +1. Los valores menores que -0.5 dejan pequeños valores negativos luego de la operación de adición. Entonces, el CTVI se realiza para resolver esta situación al dividir (NDVI + 0,50) por su valor absoluto y multiplicar el resultado por la raíz cuadrada del valor absoluto. Esto suprime el NDVI negativo. La ecuación se escribe:

Índice de vegetación transformado de Thian (TTVI)
Thiam, citado por Eastman (2003) indica que la imagen resultante del CTVI puede ser muy “ruidosa” debido a una sobrestimación de la cualidad verde. Él sugiere ignorar el primer término de la ecuación CTVI para obtener mejores resultados. Esto se logra simplemente sacando la raíz cuadrada de los valores absolutos del NDVI en la expresión original del TVI para tener un nuevo índice de vegetación llamado Índice de Vegetación Transformado de Thiam (TTVI).

Índice de vegetación de cociente (RVI)
Fue sugerido por Richardson y Wiegand, citados por  Eastman (2003) por tener gráficamente la misma fuerzas y debilidades que el TVI (ver arriba) mientras que resulta más simple a nivel computacional. El RVI es claramente el inverso del cociente simple estándar (RATIO) como se muestra en esta expresión:

Índice de vegetación cociente normalizado (NRVI)
Es una modificación del RVI por Baret y Guyot (1991) por el cual el resultado del RVI - 1 es normalizado sobre el RVI + 1.

Índice de vegetación ajustado al suelo (SAVI)
Entre los factores que modifican notablemente el comportamiento del NDVI figura la proporción de vegetación/suelo observada por el sensor. Los mismos valores de NDVI pueden corresponder a cubiertas vigorosas pero poco densas, o a cubiertas densas con poca vitalidad. Para incluir explícitamente el factor suelo, clave cuando se trabaja en zonas áridas, Huete y colaboradores propusieron incluir en la formula del NDVI un parámetro (L), que ajuste el índice a una reflectividad promedio de fondo (Huete, citado por Chuvieco 2008).

Índice de vegetación atmosféricamente resistente (ARVI)
Kaufman y Tanré, citados por Chuvieco (2008) proponen un ajuste del NDVI a las condiciones atmosféricas, teniendo en cuenta la diferente dispersión de los canales azul y rojo del espectro. De esta forma se define el ARVI de la siguiente manera:
Donde σ*IRC indica la reflectividad aparente en el infrarrojo cercano y un factor que considera la diferencia de reflectividad entre el azul y el rojo, y se define como:

Donde σ*R e indica las reflectividades aparentes en el azul y rojo, respectivamente, y ϒ es un parámetro de calibración, que depende del tipo de atmósfera, aunque para la mayor parte de los casos es igual a 1.

Índice global de monitoreo ambiental (GEMI)
Pinty y Verstraete, citados por Chuvieco (2008) proponen un índice para reducir simultáneamente el efecto atmosférico y de cambios en el color del suelo.


Índice de vegetación mejorado (EVI)
Huete, citado por Chuvieco (2008) define el índice de vegetación mejorado como una alternativa más solida a los índices tradicionales, por ser más robusto frente a la aportación del suelo y de las influencias atmosféricas. El EVI se define como:
Donde L es la corrección al efecto del fondo del follaje y C1 y C2 son coeficientes para la corrección del efecto del aerosol en las bandas rojo y azul.

Índice de infrarrojo de diferencia normalizada (NDII)
Hunt y Rock, citados por Chuvieco (2008),  proponen el NDII cuando se pretenda analizar el contenido de agua en la vegetación. Al aumentar el contenido de agua en el suelo o la vegetación, disminuye paralelamente la reflectividad en el SWIR.

Índice de área foliar (LAI)
           Esta definido por la razón entre el área foliar de toda la vegetación por unidad de área utilizada por la vegetación. El LAI es un indicador de la biomasa de cada pixel de la imagen (SEBAL 2002).

Referencias:
1.       Chuvieco Salinero, E. 2008. Teledetección Ambiental. 2 ed. Barcelona, ES, Ariel. 592 p.
2.       Eastman, R. 2003. IDRISI Kilimanjaro-Guía para SIG y Procesamiento de Imágenes. Clark Labs. 312 p.
3.       Guo, J; Mason, P. 2009. Essential Image Processing and GIS for Remote Sensing. Oxford, UK. Wiley-Blackwell.462 p.
4.       Rouse, J; Haas, R; Schell, J; Deering, D; Harlan, J. (1974). Monitoring the vernal advancement and retrogradation (Greenwave effect) of natural vegetation. Greenbelt. Maryland, US. NASA/GSFC. 87 p.
5.       SEBAL (Surface Energy Balance Algorithms for Land, US). 2002. Advanced Training and Users Manual. 97 p.
6.       Schowengerdt, R. 2007. Remote Sensing-Models and Methods for Image Processing. 3 ed. California, US, Academic Press. 558 p.
7.       Steven M. de Jong; Freek D. Van der Meer. 2004. Remote Sensing Image Analysis. California, US, Springer. 370 p.

miércoles 11 de enero de 2012

Lotka-Volterra

Modelo Presa- Depredador de Lotka-Volterra
Predator-Prey Model of Lotka Volterra
Por Eber Risco Sence
El estudio matemático de la dinámica de poblaciones data de Volterra, Lotka y Gause. Es razonable tratar el problema del modelo presa-depredador sobre las hipótesis de que el sistema, aunque muestre fluctuaciones, se mantiene en equilibrio durante cierto tiempo. Si no fuera así, el sistema ya hubiera degenerado en tiempos pasados, reduciéndose a una sola especie o a ninguna.

El modelo con ecuación diferencial más sencillo recibe el nombre de sus creadores: Lotka-Volterra. Es muy elemental, pero es un punto de partida muy útil. El modelo matemático planteado por Lotka-Volterra esta representado mediante las siguientes ecuaciones diferenciales (Pastor, 2008).
 Donde:
 : es la razón de cambio de la población de la presa con respecto al tiempo.

: es la razón de cambio de la población del depredador con respecto al tiempo.



r:  tasa per cápita de la presa.
h: probabilidad de encuentro entre el depredador y la presa y este pueda matarlo.
β: conversión eficiente de biomasa de la presa hacia la biomasa del depredador
m: probabilidad aleatoria de muerte del depredador.

Para obtener los puntos de equilibrio igualamos las dos ecuaciones  a cero. Así se tiene el nullcline para la presa (N1).
Figura 1: Nullcline para la presa (FUENTE: Pastor, 2008)

El nullcline para el depredador (N2) estaría representado de la siguiente manera:
Figura 2: Nullcline para el depredador (FUENTE: Pastor, 2008)

Además se tiene los puntos de equilibrios N1=0 y N2=0, por lo que se tiene:
Figura 3: Nullclines, equilibrios y su estabilidad, y el campo vectorial de la ecuación (FUENTE: Pastor, 2008)


Ejemplo: para r=1.1; h=0.05; β=0.2; m=0.4; población inicial de presas =60; población inicial de  predadores=40, con 10 simulaciones (población final de presas=80 y de presas=55).

  Figura 4: Modelo Lotka-Volterra para 10 simulaciones.
Figura 5: Poblaciones de las presas y depredadores con relación al tiempo.

Aplicación en MATLAB

Figura 6: Aplicación en MATLAB para dinámica poblacional.

Referencias:
1. Andrewartha, H. G.; Birch, L. C. 1954. The Distribution and Abundance of Animals. Chicago, US. University of Chicago Press.
2.Begon, M.; Harper, J. L.; Townsend, C. R. 1986. Ecology: Individuals, Populations and Communities. Oxford, UK. Blackwell.
3.  Bocarra, Nino. 2010. Modeling Complex Systems. 2 ed. New York, US. Springer. 497 p.
4. Gillma, Michael. 2009. An Introduction to Mathematical Models in Ecology and Evolution Time and Space.2 ed. Oxford, UK.167 p.
5.Lopez, José; Blé, Gamaliel. Modelo Depredador-Presa. Revista de Ciencias Básicas. UJAT (7), 25-35 (2008).
6. Murray, J.D. 2002. Mathematical Biology: An Introduction. 3 ed. New York, US. Springer. 576 p.
7. Pastor, Hohn. 2008. Mathematical Ecology of Populations and Ecosystems. Oxford, UK. Blackwell. 344 p.
8. Ranta, Esa; Lundberg, Per; kaitala, Veijo. 2006. Ecology of Populations. New York, US. Cambridge University Press. 389 p.
9. Universidad de Jaén. 2009. Modelos Matemáticos en Biología. Andalucia, ES. 360 p.

Video de Aplicación MATLAB:


Lotka Volterra por eber23

miércoles 7 de diciembre de 2011

Análisis Geomorfológico e Hidrológico de Modelos Digitales de Elevaciones con MATLAB

Análisis Geomorfológico e Hidrológico de Modelos Digitales de Elevaciones con MATLAB
Geomorphological and Hydrological Analysis of Digital Elevation Models with MATLAB
Por Eber Risco Sence


El análisis físico del terreno a partir de los valores de elevación contenidos en las celdas del Modelo Digital del Terreno, quizás por la inherente buena disposición de los mismos para el estudio mediante procedimientos y algoritmos diversos, así como por el mayor nivel de detalle y precisión respecto a las representaciones clásicas, permite obtener resultados de mucha mayor amplitud, tanto en cantidad y calidad como en el propio significado de los mismos (Olaya 2004). 

Caracterización Matemática
Evans (1972) menciona que el terreno puede ser expresado mediante una superficie cuadrática de la forma:
 
Al hacer que esta función se ajuste en la mayor medida posible a la superficie real del terreno tal y como ésta se encuentra recogida en el MDT con una menor o mayor precisión, resulta poco lógico, intuyéndose de antemano que la fidelidad para con el relieve real de tal función será poco menos que nula. No obstante, si en lugar de trabajar a nivel global con la totalidad de la malla lo hacemos a nivel local tomando el entorno inmediato de un punto, esta simplificación cobra mayor sentido y los parámetros que puedan calcularse a partir de la ecuación que se deduzca en dicha porción de la malla serán más representativos de las características reales del punto tomado (Olaya 2004).
De esa manera podemos determinar el nivel local a una malla de 9 celdas centrada alrededor de la celda denominada z5.
 
A partir de lo anterior la definición de los parámetros que configuran la superficie cuadrática a la que se pretende aproximar el entorno local de z5 de acuerdo con las siguiente expresiones, todas ellas en función de las 9 celdas consideradas en dicho entorno y de la distancia entre ellas, que coincide lógicamente con el tamaño de celda y que representamos de forma simbólica como delta de s.


De acuerdo a la formula polinómica, los coeficientes p, q, r, s ,t representan las siguientes derivadas parciales:
Figura 1: Modelo digital de elevaciones


Figura 2: Modelo digital de elevaciones mostrando curvas de nivel

Figura 3: Modelo digital de elevaciones en 3D


Parámetros geomorfológicos

Pendiente
Partiendo de la función matemática descrita anteriormente, y a la cual asimilábamos la forma del MDE en un entorno local de un punto dado, la pendiente puede calcularse a través de las primeras derivadas de dicha función (Olaya 2004).
Así la pendiente se define de la siguiente manera:
La pendiente deseada representa el ángulo entre dicho vector gradiente y la vertical, el cual, haciendo uso de conceptos básicos de cálculo de ángulos entre vectores, se obtiene según la expresión
Figura 4: Modelo digital de pendientes

Orientación
La orientación en un punto puede definirse como el ángulo existente entre el vector que señala el norte y la proyección sobre el plano horizontal del vector gradiente (Felicísimo 1994).
La forma habitual de utilizar este parámetro es expresado en grados sexagesimales, de acuerdo a la siguiente expresión (Olaya 2004).

 
Figura 5: Modelo digital de orientaciones
Sombreado
El sombreado en una superficie se obtiene de una iluminación hipotética, mediante la determinación de los valores de iluminación para cada celda. Esto se logra mediante el establecimiento de una fuente de luz hipotética y el cálculo de los valores de iluminación de cada celda en relación a las celdas vecinas (ESRI 2008).
Figura 6: Modelo digital de sombreado

Rugosidad
Upton y Fingleton, citados por Felicísimo (1994) mencionan que las coordenadas rectangulares de un vector unitario perpendicular a la superficie en un punto i vienen dadas por las expresiones:

El módulo del vector suma de un conjunto de vectores es un indicador de agrupación y, por tanto, inversamente proporcional a la rugosidad. El módulo del vector suma se calcula, para un conjunto de n datos vecinos al punto problema (los 8 más próximos, por ejemplo) mediante la expresión:
Band, citado por Felicísimo (1994) menciona que resulta conveniente estandarizar el valor de R dividiéndolo por el tamaño muestral obteniendo así el módulo medio. El resultado puede variar entre los valores extremos de 0 (dispersión máxima) y 1 (alineamiento completo). El módulo medio es complementario del parámetro estadístico denominado varianza esférica.
 Figura 7: Modelo digital de rugosidad
 Índice Elevación-Relieve
Olaya (2010) en la ayuda del programa SEXTANTE para gvSIG describe este índice según la siguiente expresión:
 Siendo Zm, Zmin, Zmax los valores medio, mínimo y máximo de elevación en la ventana de análisis 3x3 alrededor de cada celda.
 Figura 8: Índice de elevación-relieve
Índice Varianza Pendiente
Felicísimo (1994) describe este índice como el análisis mediante una ventana móvil sobre el modelo digital de pendientes. El cálculo de la varianza se realiza sobre los valores incluidos en ella y el resultado puede utilizarse como estimación de la rugosidad.
 Figura 9: Índice de varianza pendiente


Parámetros geomorfológicos (Concavidad/Convexidad)
Curvatura
Felicísimo (1994) describe la curvatura en un punto, h, puede definirse como la tasa de cambio en la pendiente y depende, por tanto, de las derivadas de segundo grado de la altitud (es decir, de los cambios de pendiente en el entorno del punto).
 Figura 10: Curvatura


Figura 11: Fuente: Hengl  y Reuter (2008)

Curvatura vertical
Krcho y Young, citados por Hengl  y Reuter (2008) definen la curvatura vertical de acuerdo a la siguiente expresión:

Figura 12: Curvatura vertical

Curvatura tangencial (horizontal)
Krcho y Hofierka, citados por Hengl  y Reuter (2008) definen la curvatura tangencial de acuerdo a la siguiente expresión:
Figura 13: Curvatura tangencial

Curvatura planta
Esta curvatura es perpendicular a la dirección de la máxima pendiente. Evans, citado por Hengl  y Reuter (2008) definen esta curvatura con la siguiente expresión:
Figura 14: Plan curvature

Curvatura inesférica (Unsphericity curvature)
Shary, citado por Hengl  y Reuter (2008) introduce un sistema de 12 curvaturas, donde las tres curvaturas básicas son: mean, unsphericity y difference curvaturas, son usadas para derivar las restantes nueve curvaturas. La Unsphericity curvatura (UNSPHC) puede definirse como:
Estos parámetros describen como la superficie difiere de una esfera y consecuentemente tiene valores de cero solo si la superficie es una esfera. 

Figura 15: Unsphericity curvature

Curvatura diferencial (Difference curvature)
Es la segunda curva independiente de las doce curvas, se define como la mitad de la diferencia entre la curvatura vertical y horizontal (Shary, citado por Hengl  y Reuter 2008).
DIFFC indica cuál de las dos curvas es formalmente la más fuerte.
  Figura 16: Difference curvature

Las tres curvas independientes (MEANC, DIFFC y UNSPHC), las restantes curvas pueden ser derivadas: la curvatura mínima y máxima (MINC, MAXC), curvatura horizontal y vertical (HEXC, VEXC), la curvatura Gaussiana total (TOTGC), curvatura acumulación total (TOTAC) y curvatura anillo total (TOTRC) (Shary, citado por Hengl  y Reuter 2008).

Curvatura Gaussiana total (Total Gaussian curvature)
Puede determinarse usando la siguiente expresión (Gauss, citado por Hengl  y Reuter 2008).
Figura 17: Total Gaussian curvature

Curvatura ROTOR (ROTOR curvature)
Florinsky, citado por Hengl  y Reuter (2008) usa también ROTOR (una curvatura de las líneas de flujo que son perpendiculares a los contornos).

ROTOR es una medida de la torsión de las líneas de flujo. Las líneas de flujo giran hacia la derecha (sentido de las agujas del reloj) cuando ROTOR>0, mientras las líneas de flujo giran en sentido anti-horario cuando ROTOR<0.
 Figura 18: ROTOR curvature

Parámetros hidrológicos
Acumulación de flujo
Hengl  y Reuter (2008) describen el principio de la acumulación de flujo: cuando la proporción d de drenaje de una celda a sus vecinos (debe sumar 1) son conocidos, también la proporción de recepción r que drena hacia una celda. Las proporciones de recepción determinan que fracción de cada vecino son recibidas. La cantidad de masa (o volumen, área u otra propiedad) A que se acumula en la celda i está dada por la suma de A en cada celda vecina multiplicada por la respectiva fracción de recepción r más la masa u otra cualidad ingresante I en la celda misma.
 Figura 19: Flujo acumulado (Matriz simple)

Figura 20: Flujo acumulado (Matriz múltiple)

Cuencas hidrográficas
La dirección de flujo es el elemento clave para la delimitación de las cuencas hidrográficas, es con el cual se da forma con toda precisión a las cuencas vertientes.
El conocimiento de esta dirección de flujo de cada celda nos permitirá deducir para cada una de ellas, y siguiendo el recorrido desde las mismas hacia aguas abajo, si pertenecen o no a la cuenca que pretendemos definir (Olaya 2004).
 Figura 21: Cuencas hidrográficas

Índice de humedad
El indice de humedad para cada celda está representado mediante la siguiente expresión (Trauth 2010)
 
 Figura 22: Índice de humedad

Índice de potencia de cauce
El índice de potencia de cauce es frecuentemente usado en geomorfología, ciencias del suelo y disciplinas relacionadas. Como una media de potencia de cauce es un indicador del potencial transporte de sedimentos y de erosión del agua. Se define mediante la siguiente expresión ((Trauth 2010)
 Figura 23: Índice de potencia de cauce

Parámetros geomorfológicos (Factores de pendiente)
Factor de profundidad S
Kinnell (s.f.) menciona que el factor de profundidad de pendiente puede ser estimado mediante la siguiente formula:
 Donde β es el ángulo de la pendiente.
Nearing, citado por Sanchez (s.f.) propone la siguiente fórmula para determinar el factor de  la profundidad de la pendiente.
 Donde θ es el ángulo de la pendiente.
Figura 24: Factor profundidad de pendiente
 Figura 25: Factor profundidad de pendiente (Nearing)
Factor Longitud L
Kinnell (s.f.) menciona que el factor de longitud de pendiente puede ser estimado mediante la siguiente formula:
Donde  Ai,j-in es el área contribuyente de escorrentía para la celda i,j y D es el tamaño  de celda, x es un factor que determina las variaciones en el ancho del flujo producto de la orientación de la celda con respecto a las curvas de nivel, m es el exponente de la longitud de la pendiente, definido por: m=0.1342*Ln(θ)+0.192
 Figura 26: Factor longitud de pendiente

Factor longitud y profundidad de Pendiente LS 
Se representa mediante el producto de los dos factores anteriores.
 
Figura 27: Factor profundidad- longitud de pendiente

Parámetros hidro-geomorfológicos (Transporte de sedimentos)
Capacidad de transporte de sedimentos
Es un índice utilizado para estimar el potencial topográfico para la erosión o deposición por medio de una expresión que representa el cambio en la capacidad de transporte de sedimentos en la dirección del flujo (Martínez 1999).
 Donde
LS=índice de capacidad de transporte de sedimentos.
AS=área de drenaje específica.
β = ángulo de la pendiente.
 Figura 28: Capacidad de transporte de sedimentos

Referencias:
1. Chavez, Richard. Sf. Modelo de Riesgo de Erosión de Suelo en cuatro micro-cuencas de la Cordillera los Maribios, Nicaragua.
2.    ESRI (Environmental Systems Research Institute, US). 2008. Manual de Ayuda de ArcGIS.
3.  Evans, I. S., General geomorphometry, derivatives of altitude, and descriptive statistics. En Chorley, R. J. 1972. (ed.) Spatial Analysis in Geomorphology, Methuen, London. pp.17-90.
4.   Felicísimo, Ángel. 1994. Modelos Digitales del Terreno: Introducción y aplicaciones en las ciencias ambientales. Pentalfa. 122 p.
5. Hengl, T., Reuter, H.I. (eds) 2008. Geomorphometry: Concepts, Software, Applications. Developments in Soil Science, vol. 33, Elsevier, 772 p.
6.    Kinnell, P. sf. The miscalculation of the USLE topographic factors in GIS. University of Camberra, Australia.
7.     Martínez-Casasnovas, J.A., 1999. Modelos digitales de terreno: Estructuras de datos y aplicaciones en el análisis de formas del terreno y en Edafología. QUADERNS DMACS Núm. 25, Departament de Medi Ambient i Ciències del Sòl, Universitat de Lleida, Lleida.
8.  Olaya Ferrero, V.2004. Hidrología Computacional y Modelos Digitales del Terreno. Madrid, ES, Víctor Olaya.365 p.
9.   Olaya Ferrero. 2010. Manual de Ayuda de SEXTANTE para gvSIG. 

Video Aplicación en MATLAB:

Análisis Geomorfológico e Hidrológico de un DEM por eber23