[GIS] How to calculate forest loss in google earth engine

change detectiongoogle-earth-engine

I am using the method outlined in a previous question – Extracting land cover type in Google Earth Engine? – to identify habitat type and area with a defined polygon. I would now like to take this information from different years to calculate forest cover loss over time. How is this best achieved?
I am using google earth engine for this process.

Best Answer

You can find global forest loss/gain from 2000-2015 in the following raster:

var hansen = ee.Image("UMD/hansen/global_forest_change_2015_v1_3");

How you can use it is explained in this Google Earth Engine tutorial.

Related Solutions

[GIS] Extracting land cover type in Google Earth Engine

Here are two ways of doing this. The first is your method of masking with some simplified code. The second is using a grouped reducer:

var s1 = ee.ImageCollection("COPERNICUS/S1_GRD"),
S2 = ee.ImageCollection("COPERNICUS/S2"),
fc = ee.FeatureCollection("ft:1w_ciwxR9FKm0rzMH8oJ3MdyE-1U15PfiOFlGW5Ib"),
fc2 = ee.FeatureCollection("ft:18Xvf7jG2EI4mJMwj0oLV2u-05bAPiwJDaxSRy2ep");

var aoi = fc; 

var bands = ['B2', 'B3', 'B4', 'B8', 'B11', 'B12'];

var S2 = ee.ImageCollection('COPERNICUS/S2')
  .filterDate('2016-07-01', '2016-08-01')
  .filterBounds(aoi);

// You may want to use median here if there's overlapping imagery.
var image = S2.mosaic();

var geometry = ee.Geometry.Rectangle([36.977296602,-15.364411987,37.057018981,-15.447436229]);
Map.centerObject(geometry);

var landcover_roi = image.clip(geometry);

Map.addLayer(landcover_roi, {'bands': 'B12,B8,B4','min': 0,'max': 4000}, 'Sentinel 2A SWIR,NIT,R', 0);

var training = landcover_roi.sampleRegions(fc2, ["info"], 30);
var classifier = ee.Classifier.cart()
    .train(training, "info", image.bandNames());
var out = landcover_roi.classify(classifier);
Map.addLayer(out, {palette: [
  'ff0000',// cleared red
  'ffcc66',// cultivated
  '387242',// forest 
  '99ff33',// grass 
  'e0e0d1',// rock 
  'e0e0d1'// rock grass
], min:1, max:6});

Map.addLayer(ee.Image().paint(fc, 1, 3));
var fc = fc2;
Map.addLayer(fc2);

var names = ['cleared', 'cultivated1', 'cultivated2', 'grass', 'rock', 'rockgrass'];
var total = out.eq([1, 2, 3, 4, 5, 6]).rename(names);

var count = total.multiply(ee.Image.pixelArea());

// The grouped reducer ignores the weights.  Unweight here for comparison.
var area = count.reduceRegion(ee.Reducer.sum().unweighted(), geometry, 30);
print(area);

var stats = ee.Image.pixelArea().addBands(out).reduceRegion({
  reducer: ee.Reducer.sum().group(1), 
  geometry: geometry, 
  scale: 30,
});
print(stats);

Google Earth Engine – Iterating Over Years for Feature Collection

You are filtering year in the correct way. This is how I'd do it:

//Load and filter the Hansen data
var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2015')
              .select(['treecover2000','loss','gain','lossyear']);

// list for filter iteration
var years = ee.List.sequence(1, 14)

// turn your scale into a var in case you want to change it
var scale = gfc2014.projection().nominalScale()

//add country districts as a feature collection
var distr = ee.FeatureCollection('ft:1U7sXFHXtxQ--g7XMeXlvPhNXPBcDtPg8Yzr2pvsg', 'geometry');

//look at tree cover, find the area
var treeCover = gfc2014.select(['treecover2000']);

// most recent version of Hansen's data has the treecover2000 layer
// ranging from 0-100. It needs to be divided by 100 if ones wants
// to calculate the areas in ha and not hundreds of ha. If not, the
// layers areaLoss/areaGain are not comparable to the areaCover. Thus
treeCover = treeCover.divide(100); // Thanks to Bruno

var areaCover = treeCover.multiply(ee.Image.pixelArea())
                .divide(10000).select([0],["areacover"])

// total loss area
var loss = gfc2014.select(['loss']);
var areaLoss = loss.gt(0).multiply(ee.Image.pixelArea()).multiply(treeCover)
               .divide(10000).select([0],["arealoss"]);

// total gain area
var gain = gfc2014.select(['gain'])
var areaGain = gain.gt(0).multiply(ee.Image.pixelArea()).multiply(treeCover)
               .divide(10000).select([0],["areagain"]);

// final image
var total = gfc2014.addBands(areaCover)
            .addBands(areaLoss)
            .addBands(areaGain)

Map.addLayer(total,{},"total")

// Map cover area per feature
var districtSums = areaCover.reduceRegions({
  collection: distr,
  reducer: ee.Reducer.sum(),
  scale: scale,
});


var addVar = function(feature) {

  // function to iterate over the sequence of years
  var addVarYear = function(year, feat) {
    // cast var
    year = ee.Number(year).toInt()
    feat = ee.Feature(feat)

    // actual year to write as property
    var actual_year = ee.Number(2000).add(year)

    // filter year:
    // 1st: get mask
    var filtered = total.select("lossyear").eq(year)
    // 2nd: apply mask
    filtered = total.updateMask(filtered)

    // reduce variables over the feature
    var reduc = filtered.reduceRegion({
      geometry: feature.geometry(),
      reducer: ee.Reducer.sum(),
      scale: scale,
      maxPixels: 1e13
    })

    // get results
    var loss = ee.Number(reduc.get("arealoss"))
    var gain = ee.Number(reduc.get("areagain"))

    // set names
    var nameloss = ee.String("loss_").cat(actual_year)
    var namegain = ee.String("gain_").cat(actual_year)

    // alternative 1: set property only if change greater than 0
    var cond = loss.gt(0).or(gain.gt(0))
    return ee.Algorithms.If(cond, 
                            feat.set(nameloss, loss, namegain, gain),
                            feat)

    // alternative 2: always set property
    // set properties to the feature
    // return feat.set(nameloss, loss, namegain, gain)
  }

  // iterate over the sequence
  var newfeat = ee.Feature(years.iterate(addVarYear, feature))

  // return feature with new properties
  return newfeat
}

// Map over the FeatureCollection
var areas = districtSums.map(addVar);

Map.addLayer(areas, {}, "areas")

In that script you get 3 fields: loss_{year}, gain_{year}, sum But if you want better 4 fields: loss, gain, year, sum; change for:

return ee.Algorithms.If(cond, 
                        feat.set("loss", loss, "gain", gain, "year", actual_year), 
                        feat)

You could also compute percentage and set it to the features.

Edit: Thank to @Bruno_Conte_Leite, who made me reconsider my answer, I have made some updates, the one suggested by Bruno and others.

Scale: I suggest to keep the original scale of Hansen data.
treeCover: most recent version of Hansen's data has the treecover2000 layer ranging from 0-100. It needs to be divided by 100 if ones wants to calculate the areas in ha and not hundreds of ha. (Bruno)
areaLoss and areaGain: Added .multiply(treeCover) otherwise the area would be of the whole pixel and not of the indicated percentage
maxPixels: I added maxPixels: 1e13 in the reduction

Best Answer

Related Solutions

[GIS] Extracting land cover type in Google Earth Engine

Google Earth Engine – Iterating Over Years for Feature Collection

Related Question