Начнем с расчета, необходимого для создания полосы, показывающей пиксели, в которых данные Хансена и др. показывают как потерю, так и усиление.
В наборе данных Хансена и др. есть полоса, пиксели которой равны 1 в случае потери и 0 в противном случае ( loss ), а также полоса, пиксели которой равны 1 в случае усиления и 0 в противном случае ( gain ). Чтобы создать полосу, где пиксели как в полосах loss , так и в полосах gain имеют значение 1, можно использовать логический метод and() для изображений. Метод and() вызывается аналогично image1.and(image2) и возвращает изображение, в котором пиксели равны 1, если и image1, и image2 равны 1, и 0 в остальных случаях:
Редактор кода (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');
Результат, увеличенный до Арканзаса со спутникового снимка, должен выглядеть примерно так, как показано на рисунке 1.
Объединив этот пример с результатом из предыдущего раздела , теперь можно воссоздать рисунок из начала урока:
Редактор кода (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');
Количественная оценка изменений лесного фонда в интересующем регионе
Теперь, когда вы более подробно ознакомились с полосами в наборе данных Хансена и др., мы можем использовать изученные концепции для вычисления статистики прироста и утраты лесов в интересующем нас регионе. Для этого нам понадобятся векторные данные (точки, линии и полигоны). Векторный набор данных представлен в Earth Engine как FeatureCollection . (Узнайте больше о коллекциях объектов и об импорте векторных данных .)
В этом разделе мы сравним общий объем потерь лесов, произошедших в Республике Конго в 2012 году, с объемом потерь лесов, произошедших на охраняемых территориях страны в то же время.
Как вы узнали из руководства по API Earth Engine , ключевым методом вычисления статистики в области изображения является reduceRegion() . ( Подробнее об уменьшении областей изображения читайте здесь .) Например, предположим, что нам нужно рассчитать количество пикселей, предположительно представляющих собой потерю лесов за исследуемый период. Для этого рассмотрим следующий код:
Редактор кода (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);
В этом примере используется редуктор ee.Reducer.sum() для суммирования значений пикселей в lossImage в пределах признака congo . Поскольку lossImage состоит из пикселей со значением 1 или 0 (соответственно, с потерями), сумма этих значений эквивалентна количеству пикселей с потерями в области.
К сожалению, запуск скрипта в текущем виде приводит к ошибке типа:
Максимальное количество пикселей в reduceRegion() по умолчанию составляет 10 миллионов. Это сообщение об ошибке означает, что Республика Конго покрывает около 383 миллионов пикселей Landsat. К счастью, reduceRegion() принимает множество параметров, один из которых ( maxPixels ) позволяет управлять количеством пикселей, используемых в вычислениях. Указание этого параметра позволяет вычислениям завершиться успешно:
Редактор кода (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);
Развернув объект, выведенный на консоль, вы увидите, что результат — 4897933 потерянных пикселей леса. Вы можете немного подкорректировать вывод в консоли, снабдив вывод метками и получив интересующий вас результат из словаря, возвращаемого функцией reduceRegion() :
Редактор кода (JavaScript)
print('pixels representing loss: ', stats.get('loss'));
Расчет площадей пикселей
Вы почти готовы ответить на вопрос, какая площадь была утрачена в Республике Конго и какая из них находилась в охраняемых зонах. Осталось перевести пиксели в фактическую площадь. Это преобразование важно, поскольку мы не всегда знаем размер пикселей, переданных в функцию reduceRegion() . Для вычисления площади в Earth Engine есть метод ee.Image.pixelArea() , который генерирует изображение, в котором значение каждого пикселя равно площади пикселя в квадратных метрах. Умножая изображение потерь на изображение площади и суммируя полученные результаты, мы получаем показатель площади:
Редактор кода (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');
Результатом стало 4 372 575 052 квадратных метров потерь за исследуемый период.
Теперь вы готовы ответить на вопрос в начале этого раздела: какая площадь лесов была утрачена в Республике Конго в 2012 году и какая ее часть находилась на охраняемых территориях?
Редактор кода (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' );
Результаты показывают, что из 348 036 295 квадратных метров леса, утраченного в Республике Конго в 2012 году, 11 880 976 из них находились на охраняемых территориях, как указано в таблице Всемирной базы данных по охраняемым территориям .
Единственные изменения между этим скриптом и предыдущим — добавление информации о защищённой территории и изменение сценария с учёта общих потерь на учёт потерь за 2012 год. Это потребовало двух изменений. Во-первых, появилось новое изображение lossIn2012 , где убыток был зафиксирован в 2012 году, где он равен 1, а в противном случае — 0. Во-вторых, поскольку название диапазона отличается ( lossyear вместо loss ), пришлось изменить название свойства в операторе печати.
В следующем разделе мы рассмотрим некоторые передовые методы расчета и построения графиков потерь лесов за каждый год, а не только за один год, как мы это делали в этом разделе.