Określanie zmian w lasach

Zacznijmy od obliczeń potrzebnych do utworzenia pasma, które pokazuje piksele, w których dane Hansen et al. wskazują zarówno utratę, jak i przyrost.

Zbiór danych Hansen et al. zawiera pasmo, którego piksele mają wartość 1 w miejscach, w których nastąpiła utrata, a w pozostałych przypadkach – 0 (loss), oraz pasmo, którego piksele mają wartość 1 w miejscach, w których nastąpił przyrost, a w pozostałych przypadkach – 0 (gain). Aby utworzyć pasmo, w którym piksele w pasmach loss i gain mają wartość 1, możesz użyć metody logicznej and() na obrazach. Metoda and() jest wywoływana jako image1.and(image2) i zwraca obraz, w którym piksele mają wartość 1, jeśli zarówno image1, jak i image2 mają wartość 1, a w pozostałych przypadkach – 0:

Edytor kodu (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');

Wynik, po powiększeniu stanu Arkansas w widoku satelitarnym, powinien wyglądać podobnie jak na rysunku 1.

Rysunek 1. Piksele z informacjami o utracie i przyroście lasów w Arkansas.

Łącząc ten przykład z wynikiem z poprzedniej sekcji, możesz teraz odtworzyć rysunek z początku samouczka:

Edytor kodu (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');

Określanie ilościowych zmian w lesie w wybranym regionie

Teraz, gdy znasz już lepiej pasma w zbiorze danych Hansen et al., możemy wykorzystać zdobytą wiedzę do obliczenia statystyk dotyczących przyrostu i utraty lasów w interesującym nas regionie. W tym celu musimy użyć danych wektorowych (punktów, linii i wielokątów). Zbiór danych wektorowych jest reprezentowany w Earth Engine jako FeatureCollection. (Dowiedz się więcej o zbiorach obiektów i o tym, jak importować dane wektorowe).

W tej sekcji porównamy całkowitą utratę lasów w Republice Konga w 2012 roku z utratą lasów na obszarach chronionych w tym kraju w tym samym czasie.

Jak dowiedzieliśmy się w samouczku dotyczącym interfejsu Earth Engine API, kluczową metodą obliczania statystyk w regionie obrazu jest reduceRegion(). (Więcej informacji o zmniejszaniu regionów obrazu) Załóżmy na przykład, że chcemy obliczyć szacunkową liczbę pikseli reprezentujących utratę lasów w okresie badania. W tym celu użyj tego kodu:

Edytor kodu (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);

W tym przykładzie użyto funkcji ee.Reducer.sum(), aby zsumować wartości pikseli w lossImage w ramach funkcji congo. Ponieważ macierz lossImage składa się z pikseli o wartości 1 lub 0 (odpowiednio dla utraty lub braku utraty), suma tych wartości jest równoważna liczbie pikseli utraty w regionie.

Niestety uruchomienie skryptu w obecnej postaci spowoduje błąd podobny do tego:

Słownik (błąd)
  Image.reduceRegion: zbyt wiele pikseli w regionie. Znaleziono 383702302, ale dozwolona jest tylko liczba 10000000.

Domyślna maksymalna liczba pikseli w reduceRegion() to 10 milionów. Ten komunikat o błędzie oznacza, że Republika Konga obejmuje około 383 mln pikseli Landsat. Na szczęście funkcja reduceRegion() przyjmuje wiele parametrów, z których jeden (maxPixels) pozwala kontrolować liczbę pikseli używanych w obliczeniach. Określenie tego parametru umożliwia wykonanie obliczeń:

Edytor kodu (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);

Rozwiń obiekt wydrukowany w konsoli i zobacz, że wynik to 4897933 pikseli utraconego lasu. Możesz nieco uporządkować wydruk w konsoli, dodając etykiety do danych wyjściowych i pobierając interesujący Cię wynik ze słownika zwróconego przez reduceRegion():

Edytor kodu (JavaScript)

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

Obliczanie obszarów pikseli

Już prawie możesz odpowiedzieć na pytanie, ile obszaru zostało utracone w Republice Konga i ile z tego obszaru znajdowało się na terenach chronionych. Pozostała część to przekształcenie pikseli na rzeczywistą powierzchnię. Ta konwersja jest ważna, ponieważ niekoniecznie znamy rozmiar pikseli wprowadzanych do reduceRegion(). Aby ułatwić obliczanie obszarów, Earth Engine ma metodę ee.Image.pixelArea(), która generuje obraz, w którym wartość każdego piksela to jego powierzchnia w metrach kwadratowych. Pomnożenie obrazu utraty przez obraz obszaru, a następnie zsumowanie wyniku daje nam miarę obszaru:

Edytor kodu (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');

Wynik to 4 372 575 052 m² utraconych w okresie badania.

Możesz teraz odpowiedzieć na pytanie z początku tej sekcji: ile obszarów leśnych zostało utraconych w Republice Konga w 2012 roku i jaka część tych obszarów znajdowała się na terenach chronionych?

Edytor kodu (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'
);

Dane wyjściowe wskazują, że z 348 036 295 m² lasów utraconych w Republice Konga w 2012 roku 11 880 976 m² znajdowało się na obszarach chronionych, co zostało przedstawione w tabeli Światowej bazy danych obszarów chronionych.

Jedyną zmianą w porównaniu z poprzednim skryptem jest dodanie informacji o obszarach chronionych i zmiana zakresu analizy z ogólnych strat na straty w 2012 roku. Wymagało to wprowadzenia 2 zmian. Po pierwsze, jest nowy obraz lossIn2012, który ma wartość 1 w miejscach, w których w 2012 r. odnotowano straty, a w pozostałych miejscach – 0. Po drugie, nazwa zespołu jest inna (lossyear zamiast loss), więc w instrukcji drukowania trzeba było zmienić nazwę właściwości.

W następnej sekcji omówimy zaawansowane metody obliczania i wykresy utraty lasów w każdym roku, a nie tylko w jednym roku, jak to zrobiliśmy w tej sekcji.