Una serie temporal de imágenes satelitales sin huecos no existe, al menos no en el sur de Chile; lo que sí existe, apenas abres tu primera descarga de Sentinel-2 para seguir el bosque que bordea el Río Calle-Calle, es un archivo lleno de parches vacíos donde debería haber verde. Esa es la pregunta que más me hacen cuando cuento que ando metida en esto de la teledetección con R: ¿de verdad se pueden rellenar esos huecos, o es puro cuento? La respuesta corta es que sí, pero no todos los huecos se rellenan igual, y hay un límite bastante claro antes de que dejes de interpolar datos y empieces a inventarlos.
¿Por qué mis series temporales de raster quedan llenas de huecos?
Lo primero que hay que entender es por qué pasa esto. Las imágenes de Sentinel-2 tienen fama de pasar cada cinco días por el mismo punto del planeta, y uno piensa que eso significa un mapa nuevo cada semana sin falta. La realidad en una región donde llueve la mitad del año es otra: de esas pasadas cada cinco días, con suerte una al mes queda limpia de nubes. Con Landsat la cosa se pone peor todavía, porque su revisita es de dieciséis días, así que si justo ese día hay neblina, te quedas sin ver el bosque por un mes entero. Cuando por fin bajas una imagen y le aplicas la máscara de nubes para que no te mienta, ahí aparecen los NA: esos huecos sin ningún valor donde el sensor simplemente no alcanzó a ver nada.
En mi caso trabajo siempre en el huso 18S, y si dejas esos huecos sueltos en una serie pensada para detectar cambios en el bosque a través del tiempo, el gráfico de evolución te queda dando saltos que no tienen nada que ver con lo que pasa de verdad allá afuera.
Rellenar huecos con approxNA (terra): qué hace y por qué la memoria se traba en series temporales largas
Después de dar varias vueltas por foros, y después de ordenar por fin lo que aprendí organizando mis scripts de R para no perder tiempo, encontré que el paquete terra (el que reemplazó al viejo raster y que cualquier guía te recomienda ahora) tiene una función llamada approxNA. La idea es simple: si tienes una imagen del primero de julio y otra del veinte, con nubes tapando los días de entremedio, el programa traza una interpolación lineal y calcula un valor razonable para lo que falta, como si dibujara una línea recta entre los dos puntos que sí tienes.
Nadie me lo explicó en una sala de clases: lo fui entendiendo de a poco, a prueba y error, revisando ejemplos de otra gente que ya había pasado por lo mismo con paquetes del ecosistema espacial de R.
La primera vez que corrí approxNA sobre una serie larga, me quedé mirando la barra de progreso avanzar a paso de tortuga durante casi media hora, hasta que apareció un mensaje que todavía me da un poco de susto: Error: vector memory exhausted. R se había rendido, no yo, aunque en el momento sentí que mi computador simplemente no daba para esto. Lo que pasaba en realidad es que estaba tratando de procesar toda la serie de una vez, sin cortarla en pedazos más chicos, y eso explota la memoria antes de terminar. Si te va a pasar lo mismo, vale la pena leer sobre cómo evitar errores de memoria al procesar rasters grandes en RStudio antes de lanzarte con una serie temporal completa, porque vas a terminar reiniciando la sesión más veces de las que te gustaría.
¿Cuándo interpolar es honesto y cuándo es pura adivinanza?
Esa interpolación lineal parte de un supuesto que no siempre es cierto: que el cambio entre una fecha y otra fue parejo, un poco cada día, como si el bosque se pintara de verde con cuentagotas. La naturaleza no funciona tan ordenada. Un temporal de varios días puede cambiarlo todo de golpe, o una helada puede quemar hojas de la noche a la mañana, y la línea recta que traza approxNA no tiene forma de saber eso.
Por eso, si rellenas de más, terminas suavizando la variabilidad real que el sensor sí alcanzó a captar, y ese suavizado puede inventar tendencias que no existen cuando después uses esa serie para mirar cambios de largo plazo. La regla que me ha servido para no pasarme es no interpolar huecos de más de un mes: si la brecha es corta, unos días o un par de semanas, la interpolación lineal es una aproximación razonable; si es más larga que eso, ya no estás rellenando, estás adivinando, y es mejor dejar el hueco marcado como lo que es, un vacío, que inventar un dato con apariencia de real.
Cómo saber si el relleno quedó bien, y no solo bonito
Para mí la prueba de que un relleno quedó bien no es que se vea liso en el gráfico, sino que las formas del terreno sigan teniendo sentido. Me acuerdo de la vez en que la pantalla por fin mostró algo que no era pura estática: el contorno del río Cruces dibujado en verde, recortado limpio contra el gris parejo alrededor, sin un parche negro cortándolo a la mitad. Esa es la prueba que uso ahora: si una forma que conozco del terreno (un río, un camino, el borde de un cerro) aparece completa y reconocible después de rellenar, probablemente el relleno tiene sentido; si en cambio el resultado se ve raro, con bordes que no calzan con nada que exista afuera, mejor desconfiar y revisar de nuevo antes de darlo por bueno.
Lo que esto no resuelve: clasificación, especies y estadísticas por polígono
Rellenar huecos en una serie no es lo mismo que recortar un raster a la forma exacta de un polígono, ni que clasificar coberturas vegetales, ni que modelar dónde podría estar presente una especie con algo como Maxent: son pasos distintos, cada uno con sus propios dolores de cabeza. Tampoco resuelve cómo extraer el valor de un raster en un punto exacto, cómo sacar estadísticas por polígono de bosque, cómo limpiar puntos de presencia antes de armar un modelo de nicho, cómo validar una clasificación con una matriz de confusión, ni cómo cambiarle la resolución a un raster remuestreándolo. Y tampoco es el lugar para meterse con la fórmula general detrás de índices como el NDVI o el NDWI, ni con el diseño de un mapa con la calidad que pide un informe técnico. Cada uno de esos temas se pelea aparte, un sábado a la vez, no todos de una.
¿Qué tengo que saber antes de llegar a este paso?
Antes de pelear con los huecos, tienes que saber leer un archivo raster en RStudio sin que se te caiga la sesión, y tienes que haber bajado tus propias imágenes de Sentinel en algún momento, ya sea desde la plataforma de la ESA o directo desde R. La primera vez que yo intenté algo parecido, bajé imágenes Landsat guardadas como archivos TIFF sin saber que esos archivos no se abren con cualquier programa; los tuve ahí guardados, sin poder verlos, durante no sé cuántos días, hasta que entendí que necesitaba algo más especializado que el visor de fotos del computador.
También ayuda tener una idea de qué son las bandas espectrales y el NDVI, porque casi siempre los huecos que vas a rellenar están en una serie de índice de vegetación y no en una imagen cualquiera, aunque ese es un tema que se explica mejor en otra parte.
Tengo una amiga que rema en kayak por el río Cruces y una vez me preguntó, medio en broma, si lo que hacía con mis mapas era dibujar el bosque o inventarlo. Le dije que ninguna de las dos cosas: lo que hago es mirar el Río Calle-Calle y sus alrededores con los huecos que el sensor dejó, y decidir, hueco por hueco, cuáles puedo rellenar con algo razonable y cuáles tengo que dejar en blanco porque ahí ni el satélite ni yo vimos lo que pasó.
Al final la regla que me sirve es bien simple: relleno los huecos cortos, los de pocos días o un par de semanas como mucho, y dejo los huecos largos tal cual están, marcados como lo que son, datos que no existen. Un hueco vacío bien marcado dice más la verdad que un parche interpolado a la fuerza, y para una serie temporal que de verdad quieras usar después, eso importa más que tener un mapa sin un solo pixel en negro.