La luz de esta tarde en Valdivia es de esas grises y parejas, sin ningún rayo directo que rebote feo contra la pantalla, perfecta para quedarse mirando un raster sin que los ojos sufran. Llevo dos años metida en esto de la teledetección con R los fines de semana, sola, sin ser programadora ni geógrafa, y el pendiente de esta vez era entender por qué el NDVI no me alcanzaba para mapear como corresponde el pedazo de bosque y jardines que hay cerca de mi casa, en el Barrio Alemania. El índice SAVI apareció como la respuesta, aunque me costó bastante entender que en RStudio hacía falta algo más que la fórmula básica de análisis espacial que ya conocía.
El NDVI no alcanzaba en un pedazo de bosque como el mío
El NDVI es lo primero que aprende cualquiera que se mete en esto, básicamente porque es la fórmula más simple y casi todos los tutoriales parten por ahí. El problema aparece cuando la vegetación no cubre completo el suelo: la tierra reflectante se cuela en el cálculo y el mapa sale con un ruido raro, como si el suelo húmedo y la maleza fueran la misma cosa. Con los renovales de mi sector, donde todavía queda harta tierra a la vista entre las plantas, el resultado se veía simplemente mal y durante varias tardes no entendí por qué.
Entendí, después de darle vueltas, que el problema no era de código sino del propio índice, había que sumarle algo más. El SAVI, o índice de vegetación ajustado al suelo, corrige justo ese ruido, aunque en su momento el nombre me sonó a otro término más de una lista larga que todavía no manejaba del todo: raster, proyección, banda espectral, cada palabra nueva era una pelea aparte antes de siquiera llegar a escribir una línea de código.
SAVI frente a NDVI: dónde está realmente la diferencia
La diferencia entre los dos no está en cuál es "mejor", sino en qué tan cerrada está la vegetación que se está mirando. Un investigador llamado Huete propuso allá por 1988 agregarle al NDVI clásico un factor de ajuste, la L, para que el suelo dejara de meter ruido en zonas donde la planta no tapa completo el terreno. Para un bosque bien cerrado, con casi nada de tierra visible, el NDVI simple ya funciona bien y complicarse con la L resulta casi innecesario; los dos índices, de hecho, salen de la misma familia de fórmulas que se calculan como una diferencia normalizada entre bandas, solo que el SAVI le suma ese ajuste extra.
Cómo armé la fórmula del SAVI en RStudio
Para calcular el SAVI hacen falta las mismas dos bandas que usa el NDVI: la del infrarrojo cercano (NIR) y la banda roja (Red). Si trabajas con imágenes de Landsat 8 o 9, que es lo que yo bajo de la página de la USGS, la NIR corresponde a la banda 5 y la Red a la banda 4, y si te equivocas de número el mapa que sale parece un cuadro abstracto que no sirve para nada.
La fórmula en sí es un poco más larga que la del NDVI clásico: ((NIR - Red) / (NIR + Red + L)) * (1 + L). Cargué la librería terra, bastante más rápida que otras opciones más antiguas para este tipo de cálculo, y definí L <- 0.5 al principio del script para no olvidarme, porque ese 0.5 es el valor estándar cuando hay una mezcla pareja de verde y tierra, como me pasa en el Barrio Alemania.
Mucho antes de manejar R medianamente bien, me acuerdo de una vez en que descargué puras imágenes Landsat como archivos TIFF pensando que las podría abrir en cualquier programa, sin saber que necesitaba algo especializado para siquiera mirarlas. Terminé esa tarde sin ver nada, con una carpeta llena de archivos que ni siquiera sabía cómo abrir.
El factor L no es un número fijo
Me quedé pensando en ese 0.5 bastante rato, porque al principio asumí que era un número fijo para siempre, como una constante universal. No es así: si estuviera mirando una selva bien cerrada, con el dosel tapando casi todo el suelo, la L debería acercarse a cero, y ahí el SAVI prácticamente se convierte en un NDVI normal. Si en cambio estuviera analizando una zona semi-desértica, con harta tierra pelada y poca planta, la L tendría que subir, más cerca de uno.
También me acuerdo de otra tarde, antes de todo esto, en que encontré sola un error de proyección que llevaba horas dándome vueltas, y de repente el mensaje de la consola, que hasta entonces me sonaba en chino, por fin tuvo sentido. Ahí entendí que el sistema de referencia de coordenadas no es un detalle menor: si la proyección no calza, ni el NDVI ni el SAVI sirven de mucho, por más bien escrita que esté la fórmula.
Cuándo quedarte con el NDVI y cuándo pasarte al SAVI
Una amiga que se la pasa remando en kayak por el río Cruces me preguntó una vez, medio en broma, para qué servía tanto verde distinto en la pantalla, y ahí caí en la cuenta de que ni yo tenía tan claro cómo explicarlo sin mirar el código de nuevo.
Si la cobertura es densa y cerrada, un bosque nativo bien tupido, por ejemplo, casi sin tierra asomando, el NDVI simple alcanza de sobra y no vale la pena complicarse con la L. Pero si el paisaje es como el mío, con parches de suelo mezclados entre la vegetación, ahí el SAVI marca la diferencia real, y conviene mover la L según qué tan cerrado o pelado esté el terreno: más cerca de cero para lo denso, más cerca de uno para lo disperso o casi desértico. Esa es, al menos, la regla que uso antes de ponerme a programar cualquier cosa.
Con el tiempo también he ido probando otras cosas sueltas los fines de semana: recortar un raster grande para quedarme solo con el pedazo que me interesa, sacarle las nubes a una imagen de Sentinel-2 antes de calcular cualquier índice, comparar escenas para ver cómo cambia el bosque de una temporada a otra, clasificar coberturas de vegetación y después validar esa clasificación con una matriz de confusión para saber si de verdad sirve, sacar el valor de un raster justo en un punto del mapa, y hasta asomarme, tímidamente, al modelado de distribución de especies con maxent. Cada una de esas cosas es un sábado entero aparte, y ni hablar de cuando el raster pesa tanto que el computador se cuelga a media tarde, ahí una aprende rápido a tener paciencia.
Después de ordenar todo esto, me dieron ganas de seguir enredándome con los pasos para realizar una interpolación espacial en R con datos climáticos, para ver si la lluvia cae igual en todo el barrio, y más adelante espero meterme con cómo crear mapas de idoneidad de hábitat en R para especies locales, aunque todavía me falta soltarme con eso de armar un mapa que se vea prolijo, con leyenda y todo, para mostrárselo a alguien más. Si estás recién empezando, mi consejo es simple: antes de copiar cualquier fórmula, fíjate primero en qué tan tapado o pelado está el suelo que estás mirando, porque de eso depende si el NDVI te alcanza o si de verdad necesitas el ajuste del SAVI.