Cuantificación del cambio forestal

Comencemos con el cálculo necesario para crear una banda que muestre los píxeles en los que los datos de Hansen et al. muestran tanto pérdida como ganancia.

El conjunto de datos de Hansen et al. tiene una banda cuyos píxeles son 1 donde se produjo pérdida y 0 en otros casos (loss), y una banda que es 1 donde se produjo ganancia y 0 en otros casos (gain). Para crear una banda en la que los píxeles de las bandas loss y gain tengan un 1, puedes usar el método lógico and() en las imágenes. El método and() se llama como image1.and(image2) y devuelve una imagen en la que los píxeles son 1 cuando tanto image1 como image2 son 1, y 0 en cualquier otro caso:

Editor de código (JavaScript)

// Load the data and select the bands of interest.
var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2015');
var lossImage = gfc2014.select(['loss']);
var gainImage = gfc2014.select(['gain']);

// Use the and() method to create the lossAndGain image.
var gainAndLoss = gainImage.and(lossImage);

// Show the loss and gain image.
Map.addLayer(gainAndLoss.updateMask(gainAndLoss),
    {palette: 'FF00FF'}, 'Gain and Loss');

El resultado, con Arkansas ampliado y la vista satelital, debería verse similar a la Figura 1.

Figura 1: Píxeles con pérdida y aumento de bosque en Arkansas.

Si combinamos este ejemplo con el resultado de la sección anterior, ahora es posible recrear la figura del comienzo del instructivo:

Editor de código (JavaScript)

// Displaying forest, loss, gain, and pixels where both loss and gain occur.
var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2015');
var lossImage = gfc2014.select(['loss']);
var gainImage = gfc2014.select(['gain']);
var treeCover = gfc2014.select(['treecover2000']);

// Use the and() method to create the lossAndGain image.
var gainAndLoss = gainImage.and(lossImage);

// Add the tree cover layer in green.
Map.addLayer(treeCover.updateMask(treeCover),
    {palette: ['000000', '00FF00'], max: 100}, 'Forest Cover');

// Add the loss layer in red.
Map.addLayer(lossImage.updateMask(lossImage),
    {palette: ['FF0000']}, 'Loss');

// Add the gain layer in blue.
Map.addLayer(gainImage.updateMask(gainImage),
    {palette: ['0000FF']}, 'Gain');

// Show the loss and gain image.
Map.addLayer(gainAndLoss.updateMask(gainAndLoss),
    {palette: 'FF00FF'}, 'Gain and Loss');

Cuantificación del cambio forestal en una región de interés

Ahora que conoces mejor las bandas del conjunto de datos de Hansen et al., podemos usar los conceptos que aprendimos hasta ahora para calcular estadísticas sobre la ganancia y la pérdida de bosques en una región de interés. Para ello, deberemos usar datos vectoriales (puntos, líneas y polígonos). Un conjunto de datos vectoriales se representa como un FeatureCollection en Earth Engine. (Obtén más información sobre las colecciones de entidades y cómo importar datos vectoriales).

En esta sección, compararemos la cantidad total de pérdida forestal que se produjo en la República del Congo en el año 2012 con la cantidad de pérdida forestal que se produjo en las áreas protegidas del país en el mismo período.

Como aprendiste en el instructivo de la API de Earth Engine, el método clave para calcular estadísticas en una región de la imagen es reduceRegion(). (Obtén más información para reducir las regiones de imágenes). Por ejemplo, supongamos que queremos calcular la cantidad de píxeles que se estima que representan la pérdida de bosques durante el período de estudio. Para ello, considera el siguiente código:

Editor de código (JavaScript)

// Load country features from Large Scale International Boundary (LSIB) dataset.
var countries = ee.FeatureCollection('USDOS/LSIB_SIMPLE/2017');

// Subset the Congo Republic feature from countries.
var congo = countries.filter(ee.Filter.eq('country_na', 'Rep of the Congo'));

// Get the forest loss image.
var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2015');
var lossImage = gfc2014.select(['loss']);

// Sum the values of forest loss pixels in the Congo Republic.
var stats = lossImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: congo,
  scale: 30
});
print(stats);

En este ejemplo, se usa el reductor ee.Reducer.sum() para sumar los valores de los píxeles en lossImage dentro del atributo congo. Dado que lossImage consta de píxeles que tienen un valor de 1 o 0 (para pérdida o no pérdida, respectivamente), la suma de estos valores equivale a la cantidad de píxeles de pérdida en la región.

Lamentablemente, ejecutar la secuencia de comandos tal como está genera un error como el siguiente:

Dictionary (Error)
  Image.reduceRegion: Demasiados píxeles en la región. Se encontró 383702302, pero solo se permiten 10000000.

La cantidad máxima predeterminada de píxeles en reduceRegion() es de 10 millones. Este mensaje de error indica que la República del Congo abarca alrededor de 383 millones de píxeles de Landsat. Por suerte, reduceRegion() toma muchos parámetros, uno de los cuales (maxPixels) te permite controlar cuántos píxeles se usan en el cálculo. Si se especifica este parámetro, el cálculo se realizará correctamente:

Editor de código (JavaScript)

// Load country features from Large Scale International Boundary (LSIB) dataset.
var countries = ee.FeatureCollection('USDOS/LSIB_SIMPLE/2017');

// Subset the Congo Republic feature from countries.
var congo = countries.filter(ee.Filter.eq('country_na', 'Rep of the Congo'));

// Get the forest loss image.
var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2015');
var lossImage = gfc2014.select(['loss']);

// Sum the values of forest loss pixels in the Congo Republic.
var stats = lossImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: congo,
  scale: 30,
  maxPixels: 1e9
});
print(stats);

Si expandes el objeto impreso en la consola, observarás que el resultado es de 4897933 píxeles de bosque perdido. Puedes limpiar un poco la impresión en la consola etiquetando el resultado y obteniendo el resultado de interés del diccionario que devuelve reduceRegion():

Editor de código (JavaScript)

print('pixels representing loss: ', stats.get('loss'));

Cómo calcular las áreas de píxeles

Ya casi puedes responder la pregunta sobre cuánta superficie se perdió en la República del Congo y cuánta de esa pérdida se produjo en áreas protegidas. La parte restante es convertir los píxeles en área real. Esta conversión es importante porque no necesariamente conocemos el tamaño de los píxeles que se ingresan en reduceRegion(). Para ayudar a calcular áreas, Earth Engine tiene el método ee.Image.pixelArea(), que genera una imagen en la que el valor de cada píxel es el área del píxel en metros cuadrados. Si multiplicamos la imagen de pérdida por esta imagen de área y, luego, sumamos el resultado, obtenemos una medida del área:

Editor de código (JavaScript)

// Load country features from Large Scale International Boundary (LSIB) dataset.
var countries = ee.FeatureCollection('USDOS/LSIB_SIMPLE/2017');

// Subset the Congo Republic feature from countries.
var congo = countries.filter(ee.Filter.eq('country_na', 'Rep of the Congo'));

// Get the forest loss image.
var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2015');
var lossImage = gfc2014.select(['loss']);
var areaImage = lossImage.multiply(ee.Image.pixelArea());

// Sum the values of forest loss pixels in the Congo Republic.
var stats = areaImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: congo,
  scale: 30,
  maxPixels: 1e9
});
print('pixels representing loss: ', stats.get('loss'), 'square meters');

Ahora, el resultado es de 4,372,575,052 metros cuadrados perdidos durante el período del estudio.

Ahora puedes responder la pregunta que se hizo al comienzo de esta sección: ¿Cuánta superficie forestal se perdió en la República del Congo en 2012 y cuánta de esa pérdida se produjo en áreas protegidas?

Editor de código (JavaScript)

// Load country features from Large Scale International Boundary (LSIB) dataset.
var countries = ee.FeatureCollection('USDOS/LSIB_SIMPLE/2017');

// Subset the Congo Republic feature from countries.
var congo = ee.Feature(
  countries
    .filter(ee.Filter.eq('country_na', 'Rep of the Congo'))
    .first()
);

// Subset protected areas to the bounds of the congo feature
// and other criteria. Clip to the intersection with congo.
var protectedAreas = ee.FeatureCollection('WCMC/WDPA/current/polygons')
  .filter(ee.Filter.and(
    ee.Filter.bounds(congo.geometry()),
    ee.Filter.neq('IUCN_CAT', 'VI'),
    ee.Filter.neq('STATUS', 'proposed'),
    ee.Filter.lt('STATUS_YR', 2010)
  ))
  .map(function(feat){
    return congo.intersection(feat);
  });

// Get the loss image.
var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2015');
var lossIn2012 = gfc2014.select(['lossyear']).eq(12);
var areaImage = lossIn2012.multiply(ee.Image.pixelArea());

// Calculate the area of loss pixels in the Congo Republic.
var stats = areaImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: congo.geometry(),
  scale: 30,
  maxPixels: 1e9
});
print(
  'Area lost in the Congo Republic:',
  stats.get('lossyear'),
  'square meters'
);

// Calculate the area of loss pixels in the protected areas.
var stats = areaImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: protectedAreas.geometry(),
  scale: 30,
  maxPixels: 1e9
});
print(
  'Area lost in protected areas:',
  stats.get('lossyear'),
  'square meters'
);

El resultado indica que, de los 348,036,295 metros cuadrados de bosque perdidos en la República del Congo en 2012, 11,880,976 se encontraban en áreas protegidas, como se representa en la tabla de la Base de datos mundial sobre áreas protegidas.

Los únicos cambios entre este y el anterior son la adición de la información del área protegida y el cambio del script para que analice la pérdida en 2012 en lugar de la pérdida general. Esto requirió dos cambios. Primero, hay una nueva imagen lossIn2012 que tiene un 1 donde se registró una pérdida en 2012 y un 0 en cualquier otro caso. En segundo lugar, debido a que el nombre de la banda es diferente (lossyear en lugar de loss), el nombre de la propiedad tuvo que cambiar en la instrucción print.

En la siguiente sección, exploraremos algunos métodos avanzados para calcular y graficar la pérdida de bosques cada año, en lugar de solo un año como lo hicimos en esta sección.