Calculating Area in Google Earth Engine

When working on Remote Sensing applications, many operations require calculating area. For example, one needs to calculate area covered by each class after supervised classification or find out how much area within a region is affected after a disaster. Calculating area for rasters and vectors is a straightforward operation in most software packages, but it is done in a slightly different way in Google Earth Engine – which can be confusing to beginners. In this post I will outline methods of calculating areas for both vectors as well as images. We will cover the following topics, starting from simple to complex.

  • Area Calculation for Features i.e. vector data
  • Area Calculation for Images (Single Class)
  • Area Calculation for Images by Class
  • Area Calculation for Images by Class by Region

Preparing the Data

Land Cover Classification

For this post, we will use the MODIS 500m Landcover dataset from the Earth Engine Data Catalog. We will select the 2018 Land Cover Image from this collection. We pick the classification scheme in the band ‘LC_Type1’ which uses 17 different land use classes. The following code snippet shows how to select the image and visualize it.

var modisLandcover = ee.ImageCollection("MODIS/006/MCD12Q1")
var filtered = modisLandcover.filter(
  ee.Filter.date('2018-01-01', '2018-12-31'))
var landcover2018 = ee.Image(filtered.first())
var classified = landcover2018.select('LC_Type1')

var palette = ['05450a', '086a10', '54a708',
 '78d203', '009900', 'c6b044','dcd159', 
 'dade48', 'fbff13', 'b6ff05', '27ff87',
 'c24f44', 'a5a5a5', 'ff6d4c', '69fff8',
 'f9ffa4', '1c0dff']
Map.addLayer(classified,
 {min:1, max:17, palette: palette},
 'MODIS Landcover 2018')

Administrative Boundaries

To compute the area for regions, we will use the FAO GAUL: Global Administrative Unit Layers 2015, Second-Level Administrative Units dataset. This is a comprehensive global dataset that contains upto Level 2 (districts/counties/…) boundaries. For this post, we will filter it and use the district boundaries for the state of Kerala in India.

var gaul = ee.FeatureCollection(
  'FAO/GAUL_SIMPLIFIED_500m/2015/level2')
var kerala = gaul.filter(
  ee.Filter.eq('ADM1_NAME', 'Kerala'))
Map.addLayer(kerala, {color: 'purple'}, 'Admin2 Boundaries')

Now that we have polygons for the region of interest, we can clip the landcover image to the extent of this collection.

var keralaLandcover = classified.clip(kerala)
Map.addLayer(keralaLandcover,
 {min:1, max:17, palette: palette},
 'Kerala Land Cover 2018')

Area Calculation for Features

Calculating area for vector features is done using the built-in area() function. A pro-tip is that you can call .geometry() on a feature collection which gives the dissolved geometry of all features in the collection. We use this method to get a geometry representing the state of Kerala and compute its area.

var stateArea = kerala.geometry().area()
var stateAreaSqKm = ee.Number(stateArea).divide(1e6).round()
print(stateAreaSqKm)
// Result: 37999

Area Calculation for Images (Single Class)

Area calculation for images is done using the ee.Image.pixelArea() function. This function creates an image where each pixel’s value is the area of the pixel. If the image pixels contains values 0 or 1 – we can multiply this pixel area image with our image and calculate the total area using reduceRegion() function.

In this example, we extract the pixels for the class Urban Areas (class 13) and see how we can calculate the total built-up area in the state.

var urban = keralaLandcover.eq(13)
Map.addLayer(urban,
 {min:0, max:1, palette: ['grey', 'blue']},
 'Built-Up')

var areaImage = urban.multiply(ee.Image.pixelArea())
var area = areaImage.reduceRegion({
  reducer: ee.Reducer.sum(),
  geometry: kerala.geometry(),
  scale: 500,
  maxPixels: 1e10
  })
var urbanAreaSqKm = ee.Number(
  area.get('LC_Type1')).divide(1e6).round()
print(urbanAreaSqKm)
// Result: 347

Area Calculation by Class

We learnt how to calculate area for a single class. But typically when you have a classified image, you want to compute area covered by each class. We must follow a similar process as before, but using a Grouped Reducer.

var areaImage = ee.Image.pixelArea().addBands(
      keralaLandcover)

var areas = areaImage.reduceRegion({
      reducer: ee.Reducer.sum().group({
      groupField: 1,
      groupName: 'class',
    }),
    geometry: kerala.geometry(),
    scale: 500,
    maxPixels: 1e10
    }); 

print(areas)

The result of reduceRegion() with a grouped reducer is a dictionary of dictionaries for each class. The top level dictionary has a single key named ‘groups‘. To extract the individual dictionaries and get properly formatted results, we must do a little post-processing. But before we dive into that, you must learn about another important function in Earth Engine called flatten().

flatten() is an important function in Earth Engine required for data processing. It takes a nested list and converts it to a flat list. Many Earth Engine constructors such a ee.Dictionary, ee.FeatureCollection etc. expect a flat list. So before creating such objects with nested objects, we must convert them to flat structures.

var nestedList = ee.List(
  [['a', 'b'], ['c', 'd'], ['e', 'f']])
print(nestedList) 
// Output: [["a","b"],["c","d"],["e","f"]]
print(nestedList.flatten())
// Output: ["a","b","c","d","e","f"]

Now we can take the results of the grouped reducer and map a function over it to extract the individual class areas and convert it to a single dictionary. Important to note that dictionary key must be of type ‘string’. Our keys are class numbers, so we use the format() method to convert the number to string

var classAreas = ee.List(areas.get('groups'))

var classAreaLists = classAreas.map(function(item) {
  var areaDict = ee.Dictionary(item)
  var classNumber = ee.Number(areaDict.get('class')).format()
  var area = ee.Number(
    areaDict.get('sum')).divide(1e6).round()
  return ee.List([classNumber, area])
})

var result = ee.Dictionary(classAreaLists.flatten())
print(result)

Area Calculation by Class by Admin Area

We saw how we can calculate areas by class for the whole state. What if we wanted to know the breakup of these classes by each district? This requires one more level of processing. We can apply a similar computation as above, but wrap the computation in a function and by apply it using map() on the Feature Collection to obtain the values by each district geometry.

var calculateClassArea = function(feature) {
    var areas = ee.Image.pixelArea().addBands(classified)
    .reduceRegion({
      reducer: ee.Reducer.sum().group({
      groupField: 1,
      groupName: 'class',
    }),
    geometry: feature.geometry(),
    scale: 500,
    maxPixels: 1e10
    })

    var classAreas = ee.List(areas.get('groups'))
    var classAreaLists = classAreas.map(function(item) {
      var areaDict = ee.Dictionary(item)
      var classNumber = ee.Number(
        areaDict.get('class')).format()
      var area = ee.Number(
        areaDict.get('sum')).divide(1e6).round()
      return ee.List([classNumber, area])
    })

    var result = ee.Dictionary(classAreaLists.flatten())
    var district = feature.get('ADM2_NAME')
    return ee.Feature(
      feature.geometry(),
      result.set('district', district))
}

var districtAreas = kerala.map(calculateClassArea);

One thing to note is that each district may or may not have all of the 17 classes present. So each feature will have different number of attributes depending on which classes are present. We can explicitly set the expected fields in the output using the selectors argument for Export.table.toDrive() function. Because we need to use the list of output fields in the Export function we have to call getInfo() to get the list values on the client-side.

This result will be homogeneous table with all classes. Once done, we can export the results to a CSV file.

var classes = ee.List.sequence(1, 17)
var outputFields = ee.List(
    ['district']).cat(classes).getInfo()

Export.table.toDrive({
    collection: districtAreas,
    description: 'class_area_by_district',
    folder: 'earthengine',
    fileNamePrefix: 'class_area_by_district',
    fileFormat: 'CSV',
    selectors: outputFields
    })

You can see the full script with comments at https://code.earthengine.google.co.in/1cff951ce00156bacb5eba920d51f3bf

Here is the resulting CSV file.

Hope you found the post useful and got some inspiration to apply it to your own area calculation problem. If you are new to Earth Engine and want to master it, check out my course Applied Remote Sensing with Google Earth Engine.

5 Comments

Leave a Comment

    • To share your code, click the ‘Get Link’ button and paste the link. Also make any assets used in the code public so it can be accessed.

Leave a Reply