Lab 24 - Watching the Forest Disappear

Before We Start

• Prerequisites: Reducers module, Time Series module, Deforestation & Forest Monitoring module.
• Estimated time: 45 minutes
• You will need: GEE Code Editor access

Let's Get Started

Part 1: What did the forest look like in 2000?

Let's start by loading the Hansen dataset and looking at the tree cover layer for the year 2000. This baseline shows us what the forest looked like before the study period began. Think of it as the "before" photo.

// ===== Part 1: Load Hansen GFC and map tree cover =====
var roi = ee.Geometry.Rectangle([-64, -12, -60, -8]); // Rondonia region

// Load the Hansen Global Forest Change dataset
var gfc = ee.Image('UMD/hansen/global_forest_change_2023_v1_11');

// Visualize tree cover in 2000
var treeCover = gfc.select('treecover2000');
Map.centerObject(roi, 8);
Map.addLayer(treeCover.clip(roi), {
  min: 0, max: 100,
  palette: ['FFFFCC', '41AB5D', '006837']
}, 'Tree Cover 2000 (%)');

// Create a binary forest mask at 50% threshold
var forest2000 = treeCover.gte(50);
Map.addLayer(forest2000.selfMask().clip(roi), {
  palette: ['006837']
}, 'Forest > 50% Canopy');

print('GFC bands:', gfc.bandNames());

Understanding the code:

• ee.Geometry.Rectangle([-64, -12, -60, -8]) defines a bounding box over Rondonia.
• treecover2000 shows canopy cover percentage at the start of the study period.
• The 50% threshold converts continuous tree cover into a binary forest/non-forest map.

Part 2: Where did the forest go?

Now let's overlay the forest loss layer on top of our tree cover map. The loss band is binary: 1 where loss occurred, 0 where it did not. Think of it like a red marker highlighting every spot where the forest was cut down.

// ===== Part 2: Map forest loss =====
var loss = gfc.select('loss');

// Show all loss as a red overlay
Map.addLayer(loss.selfMask().clip(roi), {
  palette: ['FF0000']
}, 'All Forest Loss (2001-2023)');

// Show loss only within areas that were forest in 2000
var forestLoss = loss.and(forest2000);
Map.addLayer(forestLoss.selfMask().clip(roi), {
  palette: ['CC0000']
}, 'Forest Loss (>50% baseline)');

// Also show forest gain for context
var gain = gfc.select('gain');
Map.addLayer(gain.selfMask().clip(roi), {
  palette: ['00FF00']
}, 'Forest Gain (2000-2012)');

Understanding the code:

• .and(forest2000) restricts loss to pixels that were actually forest in 2000.
• This prevents counting loss in areas that had sparse tree cover (e.g., savanna).
• The gain band only covers 2000 to 2012 and is shown for reference.

Part 3: When did it happen?

Here is where it gets really interesting. The lossyear band tells us the exact year each pixel was deforested. Mapping it with a color ramp lets us watch deforestation march across the landscape over two decades.

// ===== Part 3: Visualize loss year with color ramp =====
var lossYear = gfc.select('lossyear');

Map.addLayer(lossYear.updateMask(forestLoss).clip(roi), {
  min: 1, max: 23,
  palette: [
    '#FFFFB2', '#FED976', '#FEB24C', '#FD8D3C', '#FC4E2A',
    '#E31A1C', '#BD0026', '#800026', '#4A0048', '#0000FF'
  ]
}, 'Loss Year (2001-2023)');

// Isolate loss from the most recent 5 years (2019-2023)
var recentLoss = lossYear.gte(19).and(lossYear.lte(23)).and(forest2000);
Map.addLayer(recentLoss.selfMask().clip(roi), {
  palette: ['FF00FF']
}, 'Recent Loss (2019-2023)');

Understanding the code:

• Loss year values are years since 2000: 1 = 2001, 23 = 2023.
• The color ramp runs from yellow (early loss) to blue (recent loss).
• .gte(19).and(.lte(23)) isolates the five most recent years.

Part 4: Putting a number on the destruction

Let's convert those loss pixels into a real measurement. How many hectares of forest have we lost? We will use ee.Image.pixelArea() to find out.

// ===== Part 4: Calculate total loss area =====
var lossAreaImage = forestLoss.multiply(ee.Image.pixelArea());

var totalLoss = lossAreaImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: roi,
  scale: 30,
  maxPixels: 1e9
});

// Convert to hectares and km²
var lossHa = ee.Number(totalLoss.get('loss')).divide(1e4);
var lossKm2 = lossHa.divide(100);

print('===== Total Forest Loss =====');
print('Hectares:', lossHa);
print('Square kilometers:', lossKm2);

// Also calculate original forest area for context
var forestAreaImage = forest2000.multiply(ee.Image.pixelArea());
var totalForest = forestAreaImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: roi,
  scale: 30,
  maxPixels: 1e9
});
var forestHa = ee.Number(totalForest.get('treecover2000')).divide(1e4);
print('Original forest area (hectares):', forestHa);
print('Percentage of forest lost:', lossHa.divide(forestHa).multiply(100));

Part 5: The annual story

Now let's break our total loss down by year. We will use a grouped reducer to get the area lost each year, then chart it. This is where patterns start to jump out: which years were the worst?

// ===== Part 5: Annual deforestation chart =====
// Create image with area and loss year bands
var areaByYear = ee.Image.pixelArea()
  .addBands(lossYear)
  .updateMask(forestLoss);

// Grouped reduction: sum area for each loss year
var annualStats = areaByYear.reduceRegion({
  reducer: ee.Reducer.sum().group({
    groupField: 1,
    groupName: 'year'
  }),
  geometry: roi,
  scale: 30,
  maxPixels: 1e9
});

// Convert to feature collection for charting
var groups = ee.List(annualStats.get('groups'));
var annualData = groups.map(function(item) {
  var dict = ee.Dictionary(item);
  var year = ee.Number(dict.get('year')).add(2000);
  var areaHa = ee.Number(dict.get('sum')).divide(1e4);
  return ee.Feature(null, {year: year, area_ha: areaHa});
});
var annualFC = ee.FeatureCollection(annualData);

// Create bar chart
var chart = ui.Chart.feature.byFeature(annualFC, 'year', 'area_ha')
  .setChartType('ColumnChart')
  .setOptions({
    title: 'Annual Forest Loss in Rondonia (2001-2023)',
    hAxis: {title: 'Year', format: '####'},
    vAxis: {title: 'Forest Loss (hectares)'},
    colors: ['#CC0000'],
    legend: {position: 'none'}
  });
print(chart);

Understanding the code:

• .group({groupField: 1}) groups the sum by the second band (lossyear, index 1).
• Each group contains a year code and a summed area value.
• .add(2000) converts year codes (1 to 23) to calendar years (2001 to 2023).
• ui.Chart.feature.byFeature() creates a chart from the feature collection.

Part 6: Do roads drive deforestation?

Let's compare deforestation between two sub-regions: one along a major highway (where we expect heavy clearing) and one deeper in the interior. The contrast should tell us a lot about what drives deforestation.

// ===== Part 6: Compare two sub-regions =====
// Region A: Along the highway (active deforestation frontier)
var regionA = ee.Geometry.Rectangle([-63, -11, -61.5, -9.5]);

// Region B: Further interior (deeper forest)
var regionB = ee.Geometry.Rectangle([-61.5, -11, -60, -9.5]);

// Calculate loss for Region A
var lossA = forestLoss.multiply(ee.Image.pixelArea()).reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: regionA,
  scale: 30,
  maxPixels: 1e9
});
var lossHaA = ee.Number(lossA.get('loss')).divide(1e4);

// Calculate loss for Region B
var lossB = forestLoss.multiply(ee.Image.pixelArea()).reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: regionB,
  scale: 30,
  maxPixels: 1e9
});
var lossHaB = ee.Number(lossB.get('loss')).divide(1e4);

print('===== Sub-Region Comparison =====');
print('Region A - Highway corridor (hectares):', lossHaA);
print('Region B - Interior forest (hectares):', lossHaB);

// Visualize the two regions
Map.addLayer(ee.Image().paint(regionA, 0, 3), {palette: ['yellow']}, 'Region A (Highway)');
Map.addLayer(ee.Image().paint(regionB, 0, 3), {palette: ['cyan']}, 'Region B (Interior)');