Calcular 23 índices de vegetación usados en agricultura con imagenes Landsat y Sentinel
En esta entrada se va explicar como calcular rápidamente cerca de 23 índices de vegetación (IV) usados comúnmente en aplicaciones agrícolas ya sea para estimar variables como cobertura, índice de área foliar, biomasa, entre otros. Para ello vamos a usar la herramienta VICAL (Jiménez-Jiménez et al, 2022) y Google Earth Engine (GEE). La sección de ayuda de VICAL se encuentra aquí.
VICAL cuenta con tres archivos que puede importar a su Script, estos son:
// archivo para las colecciones de imagenes);
var imp = require('users/InifapCenidRaspa/VICAL:Exportaciones');
// archivo para los indices de vegetacion
var imp2= require('users/InifapCenidRaspa/VICAL:VegetationIndex');
// archivo para las visualizaciones
var St= require('users/InifapCenidRaspa/VICAL:Style');
Variables generales
Antes de importar el conjunto de colecciones de imágenes (Landsat o sentinel-2 o ambos) se deben definir ciertas variables que son útiles para filtrar esta colección, estos son: i) un punto o polígono; ii) *intervalo de fechas y iii) valor umbral de nubes en las imágenes*. Para este ejemplo usamos polígonos de unas parcelas agrícolas de un Distrito de Riego en México. co
var fecha = ['2021-01-01', '2022-03-18']; //Fecha inicial y final
//polígono o punto
var table = ee.FeatureCollection("projects/calcium-verbena-328905/assets/Bate");
var p_nubes= 30;//Porcentaje de nubes
Selección de Colección de imagenes
Con estas colecciones de imágenes se pueden calcular series de tiempo de diferentes índices de vegetación.
Landsat
Si desea usar imágenes LandSat (4, 5, 7, 8 y 9) corregidas atmosféricamente libre de nubes se puede usar el siguiente código, donde se crea una función para unir las colecciones de imágenes y ordenarlos por fecha. Para ello se usa el archivo imp descrito al inicio de la entrada. Se realiza el ajuste espectral de los sensores del Landsat, 4, 5 y 7.
function ColeccionImagenSR(fecha, recorte, umbral)
{
// se importan las colecciones de imagenes usando el archivo "imp"
var L9sr = imp.ColeccionLandsatSR(fecha, 'LC09', recorte, umbral);
var L8sr = imp.ColeccionLandsatSR(fecha, 'LC08', recorte, umbral);
var L7sr = imp.ColeccionLandsatSR(fecha, 'LE07', recorte, umbral);
var L5sr = imp.ColeccionLandsatSR(fecha, 'LT05', recorte, umbral);
var L4sr = imp.ColeccionLandsatSR(fecha, 'LT04', recorte, umbral);
//Los datos de ETM y ETM+ se ajustan espectrales a OLI y OLI-2
var L7a = L7sr.map(imp.TMaOLI);
var L5a = L5sr.map(imp.TMaOLI);
var L4a = L4sr.map(imp.TMaOLI);
// Une las tres series de imágenes
var serieT =L9sr.merge(L8sr).merge(L7a).merge(L5a).merge(L4a).sort('system:time_start');
return serieT;
}
//Se importa la colección usando la función anterior
var l8Sergio=ColeccionImagenSR(fecha, table, p_nubes);
//podemos imprimir las imagenes usando la funcion print() para ver
//si se llevo a cabo el filtrado de la colección de imagenes (Figura 6.1)
print (l8Sergio);
Sentinel-2
Si se desea usar imágenes Sentinel-2 corregidas atmosféricamente libre de nubes se puede usar el siguiente código.
//Se importa la colección usando el siguiente código
var S2sr = imp.ColeccionImagenSentinelSR(fecha, table, p_nubes);
//podemos imprimir las imagenes usando la funcion print() para ver
//si se llevo a cabo el filtrado de la colección de imagenes (Figura 6.2)
print (S2sr);
Landsat y Sentinel-2
Si se desea usar imágenes LandSat y Sentinel-2 corregidas atmosféricamente libre de nubes se puede usar el siguiente código, los datos se ajustaron espectralmente a las bandas espectrales de Landsat 8. Se emplean las funciones descritas en la Sección anterior:
function ColeccionImagenAMBOS(fecha, recorte, umbral)
{
//Se lee la función para imagenes Landsat con ajuste espectral
var L8Conjunto=ColeccionImagenSR(fecha, recorte, umbral)
//colección entinel
var S2sr = imp.ColeccionImagenSentinelSR(fecha, recorte, umbral);
var S2a = S2sr.map(imp.MSIaOLI); //Se hace el ajuste espectral de sentinel-2 a Landsat
var serieT = S2a.merge(L8Conjunto).sort('system:time_start');
return serieT;
}
//Se importa la colección usando el siguiente código
var S2B = ColeccionImagenAMBOS(fecha, table, p_nubes);
//podemos imprimir las imagenes usando la funcion print() para ver
//si se llevo a cabo el filtrado de la coleccion de imagenes (Figura 6.3)
print (S2sr);
Índices de vegetación
Para usar el código de índices de vegetación de VICAL se tienen que usar el archivo imp2; y se mandan a llamar estos IV usando los nombres de la columna ExpresionGEE que se muestran en la Tabla 1.
Por ejemplo, para calcular el NDVI (Fila 12 de Tabla 1)con la colección de imágenes LandSat y Sentinel-2 de la sección anterior se usaría el siguiente código:
//Normalized Difference Vegetation Index- NDVI
var ivs = ee.ImageCollection(S2B.map(imp2.NDVI));
//para imprimir y ver la banda de NDVI
print (ivs);
Ahora, para visualizar en el mapa el NDVI de la primera imagen de la colección y recortado para la zona quedaría como se muestra a continuación; donde para ver la paleta de colores se utilizaría el archivo st de VICAL.
//Se obtiene el NDVI de la primera imagen en la colección
var iv = ivs.first();
//Paleta de colores donde se usa el archivo st
var ivVis = {min :0, max : 1, palette : St.paletaIV};
Map.addLayer(iv.clip(table), ivVis,'NDVI'); //Indice
//se centra el mapa a la zona
Map.centerObject(table, 13);
En la siguiente Figura se muestra el mapa de NDVI para la zona de interés.
Para revisar el código del ejemplo clic aquí
Si se desea visualizar el NDVI de una imagen en particular se debe convertir a lista de después mandar a llamarlo.
Table 1. Nombre de los Índices de vegetación en VICAL
Número | Indice | Abreviatura | ExpresionGEE | Coeficientes |
---|---|---|---|---|
1 | Atmospherically resistant vegetation index | ARVI | ARVI | γ=1.0 |
1 | Atmospherically resistant vegetation index | ARVI | ARVI | γ=1.0 |
2 | Adjusted transformed soil-adjusted vegetation index | ATSAVI | ATSAVI | |
3 | Difference vegetation index | DVI | DVI | |
4 | Enhanced vegetation index | EVI | EVI | C1=6.0, C2= 7.5; L=1.0 |
5 | Enhanced vegetation index | EVI2 | EVI2 | C1=2.4 |
6 | Green normalized difference vegetation index | GNDVI | GNDVI | |
7 | Modified soil adjusted vegetation index | MSAVI2 | MSAVI2 | |
8 | Moisture stress index | MSI | MSI | |
9 | Modified triangular vegetation index | MTVI | MTVI | |
10 | Modified triangular vegetation index-2 | MTVI2 | MTVI2 | |
11 | Normalized difference tillage index (NDTI) | NDTI | NDTI | |
12 | Normalized difference vegetation index | NDVI | NDVI | |
13 | Normalized difference water index | NDWI | NDWI | |
14 | Optimized soil adjusted vegetation index | OSAVI | OSAVI | X=0.16 |
15 | Renormalized difference vegetation index | RDVI | RDVI | |
16 | Redness index | RI | RI | |
17 | Ratio vegetation index | RVI | RVI | |
18 | Soil adjusted vegetation index | SAVI | SAVI | L=0.5 |
19 | Triangular vegetation index | TVI | TVI | |
20 | Transformed soil adjusted vegetation index | TSAVI | TSAVI | a= 1 ; b=0; |
21 | Visible atmospherically resistant index | VARI | VARI | |
22 | Vegetation index number or simple ratio | VIN | VIN | |
23 | Wide dynamic range vegetation index | WDRVI | WDRVI | α=0.2 |