Proyecto Peridotita 7: Poniendo a prueba nuestra clasificación — comparándola con los datos oficiales de Copernicus

Introducción

Hasta ahora, en el Proyecto Peridotita hemos usado nuestros propios umbrales de NBR para clasificar la severidad del incendio — unos umbrales razonables, pero inventados por nosotros. En esta entrada damos un paso más serio: en vez de comparar nuestra clasificación casera contra un mapa oficial hecho con otro método, reproducimos la metodología estándar internacional (el dNBR, con los umbrales del Servicio Geológico de Estados Unidos, USGS) y comprobamos si nuestro resultado coincide con el informe oficial que la Unión Europea elaboró sobre este mismo incendio a través de su servicio de emergencias Copernicus EMS.

🔬 Pregunta de investigación

Aplicando la misma metodología estándar (dNBR) que usan las agencias de emergencia profesionales, ¿cuántas hectáreas de Sierra Bermeja clasificamos como quemadas, y qué tan cerca queda esa cifra de la evaluación oficial de Copernicus?

Requisitos previos

• Vuestro polígono de área de estudio ya dibujado en Google Earth Engine.
• Una cuenta de Google Earth Engine con proyecto activo (para poder subir assets).

Paso 1: La metodología estándar — dNBR

En vez del NBR de una sola fecha, la comunidad científica usa el dNBR (NBR diferencial): se resta el NBR de antes del incendio menos el de después. Cuanto más alto y positivo el resultado, más severo el quemado. El USGS estableció hace años una tabla de umbrales que se ha convertido en el estándar de facto:

dNBR	Categoría
< -0,25	Regeneración vigorosa
-0,25 a -0,10	Regeneración leve
-0,10 a +0,10	Sin quemar
+0,10 a +0,27	Severidad baja
+0,27 a +0,44	Severidad moderada-baja
+0,44 a +0,66	Severidad moderada-alta
> +0,66	Severidad alta

Paso 2: Elegir las fechas correctas

Para que la comparación con Copernicus sea justa, hay que usar imágenes de fechas parecidas a las que ellos usaron. Según su propio informe, evaluaron la situación con una imagen del 18 de septiembre de 2021 (10 días después de que empezara el incendio, el 8 de septiembre). Por eso usamos:

• Antes del incendio: 29 de agosto de 2021
• Después del incendio: 18 de septiembre de 2021 (en vez de una fecha más tardía, para que el rebrote inicial de matorral no suavice artificialmente la severidad medida)

Paso 3: Calcular el dNBR en Google Earth Engine

(Pega esto en el editor de código de Earth Engine, con vuestro polígono geometry ya dibujado)

// ================================================================
// PROYECTO PERIDOTITA - dNBR ESTÁNDAR (metodología USGS)
// ================================================================

if (typeof geometry === 'undefined') {
  print('❌ ERROR: dibuja el polígono en el mapa antes de darle a Run.');
} else {

var aoi = geometry;
Map.centerObject(aoi, 12);

var fechaPre = ee.Date('2021-08-29');
var fechaPost = ee.Date('2021-09-18');

function imagenMasLimpia(fecha) {
  return ee.ImageCollection('COPERNICUS/S2_SR_HARMONIZED')
    .filterBounds(aoi)
    .filterDate(fecha, fecha.advance(10, 'day'))
    .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 20))
    .sort('CLOUDY_PIXEL_PERCENTAGE')
    .first();
}

var imgPre = ee.Image(imagenMasLimpia(fechaPre));
var imgPost = ee.Image(imagenMasLimpia(fechaPost));

var nbrPre = imgPre.normalizedDifference(['B8', 'B12']).rename('NBR_pre');
var nbrPost = imgPost.normalizedDifference(['B8', 'B12']).rename('NBR_post');
var dNBR = nbrPre.subtract(nbrPost).rename('dNBR');

var clasificado = dNBR.expression(
  "b('dNBR') < -0.25 ? 0" +
  " : b('dNBR') < -0.10 ? 1" +
  " : b('dNBR') < 0.10 ? 2" +
  " : b('dNBR') < 0.27 ? 3" +
  " : b('dNBR') < 0.44 ? 4" +
  " : b('dNBR') < 0.66 ? 5" +
  " : 6"
);

var paletaUSGS = ['7A8737', 'ACBE4D', 'FFFFFF', 'FFF7BC', 'FEC44F', 'D95F0E', '7F0000'];

Map.addLayer(clasificado.clip(aoi), {min: 0, max: 6, palette: paletaUSGS},
  'dNBR clasificado (metodología USGS)');

var nombresClases = ['Regeneración vigorosa', 'Regeneración leve', 'Sin quemar',
  'Severidad baja', 'Severidad moderada-baja', 'Severidad moderada-alta', 'Severidad alta'];

var areaImagen = ee.Image.pixelArea().divide(10000);

for (var i = 0; i <= 6; i++) {
  var mascara = clasificado.eq(i);
  var area = areaImagen.updateMask(mascara).reduceRegion({
    reducer: ee.Reducer.sum(), geometry: aoi, scale: 20, maxPixels: 1e9
  });
  print(nombresClases[i] + ':', area.get('area'), 'hectáreas');
}

}

Para verificar este paso: deberíais ver en el mapa un patrón muy parecido al de vuestros mapas de NBR anteriores — rojo/naranja en el centro de la mancha quemada, verde en los bordes de regeneración temprana.

Paso 4: Conseguir el mapa oficial de Copernicus

1. Entra en emergency.copernicus.eu/EMSR545 y busca el “Grading Product” (producto de clasificación de severidad).
2. Descarga el paquete vectorial (.zip) — dentro hay varios grupos de archivos: AreaOfInterestA, BuiltUpA, HydrographyA, ImageFootprintA, NaturalLandUseA, ObservedEventA, TransportationL… El que necesitáis es NaturalLandUseA — es el que clasifica la vegetación por severidad de daño, igual que hacéis vosotros con el dNBR.
3. Descomprime y quédate con los 4 archivos de NaturalLandUseA: .shp, .shx, .dbf, .prj.

Paso 5: Subir el shapefile a Earth Engine como “asset”

1. En code.earthengine.google.com, ve a la pestaña “Assets” (panel izquierdo) → “NEW” → “Shape files (.shp, .shx, .dbf, .prj)”.
2. Sube los 4 archivos juntos, ponle un nombre (por ejemplo NaturalLandUseA), y espera a que termine de procesarse.
3. Copia la ruta del asset (los tres puntos junto al nombre → copiar ruta) — la necesitaréis en el siguiente script.

Para verificar este paso: el asset debe aparecer con un icono verde en la lista, no un reloj de espera ni un aviso rojo.

Paso 6: Sacar las hectáreas oficiales por categoría

var perimetroOficial = ee.FeatureCollection('RUTA_DE_TU_ASSET');

// Primero, inspecciona los nombres de columna reales
print('Primer elemento:', perimetroOficial.first());

// El campo de severidad se llama "damage_gra", con valores:
// "Destroyed", "Damaged", "Possibly damaged"
var categorias = ['Destroyed', 'Damaged', 'Possibly damaged'];
var colores = ['7F0000', 'D95F0E', 'FFF7BC'];

categorias.forEach(function (cat, i) {
  var subset = perimetroOficial.filter(ee.Filter.eq('damage_gra', cat));
  Map.addLayer(subset, {color: colores[i]}, 'Oficial: ' + cat);
  var areaHa = subset.geometry().area().divide(10000);
  print('Hectáreas oficiales en "' + cat + '":', areaHa);
});

var quemadoTotal = perimetroOficial.filter(ee.Filter.inList('damage_gra', categorias));
var areaTotalHa = quemadoTotal.geometry().area().divide(10000);
print('🔥 TOTAL quemado (Destroyed+Damaged+Possibly damaged):', areaTotalHa, 'ha');

Para verificar este paso: el nombre del campo puede variar ligeramente según la versión del producto — si damage_gra no existe, mira el resultado de “Primer elemento” y busca la columna cuyos valores sean “Destroyed”/”Damaged”/”Possibly damaged”.

⚠️ Aviso práctico: si intentáis exportar esta capa oficial como imagen con Export.image.toDrive(), es fácil toparse con el error Request payload size exceeds the limit: 10485760 bytes — el shapefile oficial tiene miles de polígonos muy detallados (en nuestro caso, 6.944), y convertirlo directamente a imagen sobrecarga la petición. La solución más simple y rápida: una captura de pantalla del mapa ya coloreado en Earth Engine sirve perfectamente para comparar visualmente en el blog, sin necesidad de resolver ese error.

Resultados

Con nuestros datos reales de Sierra Bermeja:

Categoría oficial (Copernicus)	Hectáreas
Destroyed	102,3
Damaged	5.780,4
Possibly damaged	2.793,1
Total oficial	8.675,8 ha
Total según el informe del PDF	8.679,3 ha

La diferencia entre nuestra suma calculada en Earth Engine y la cifra que da el propio informe de Copernicus es de solo 3,5 hectáreas — un 0,04%, prácticamente idéntica. Esto confirma que estamos leyendo e interpretando correctamente los datos oficiales.

El siguiente paso natural — comparar esta cifra oficial con la que sale de nuestra propia clasificación por dNBR (Paso 3), en vez de solo con los datos en bruto de Copernicus — es donde de verdad se pone a prueba si nuestra metodología casera se parece a la profesional. (Aquí se completa con vuestra propia cifra de hectáreas quemadas del Paso 3, sumando las categorías de severidad baja, moderada-baja, moderada-alta y alta.)

Conclusión

Reproducir exactamente el mismo dato que ya tiene Copernicus (los 8.675,8 ha frente a las 8.679,3 oficiales) demuestra que sabemos leer e interpretar correctamente productos de datos científicos reales — un paso importante antes de fiarnos de nuestra propia clasificación por dNBR para el resto del proyecto. La comparación fina entre “cuánto coincide nuestro propio método con el oficial” (no solo si sabemos leer sus datos) la dejamos para trabajarla en clase con todo el grupo, junto con la discusión sobre qué significan las diferencias que puedan aparecer entre categorías.

---

Etiquetas: AIDARAC, Climate Detectives, Ciencia de Datos, Copernicus EMS, dNBR, ESA, EMSR545, ies monterroso, Málaga, montesteam, nbr, Proyecto Peridotita, sierra bermeja, Teledetección, USGS