2

How do you calculate the number of occurrences of a category that occur in each cluster using Google Earth Engine? My goal is to figure out the percentage (or # of pixels) of each landcover type that occur within each patch.

The script below shows a working example using a modified version of the developers guide (https://developers.google.com/earth-engine/image_objects and https://developers.google.com/earth-engine/reducers_grouping) to create clusters with an unique id and its corresponding landcover classification. However when I group a count reducer on the image the patch id is not correct.

The ultimate goal is to figure out at what point if one adds the percentage of landcover types sequentially (e.g. landcover 0 + 1 + 2...) does the proportion exceed 50% (i.e. the proportion >50% at landcover 9 for patch #19). Is there a better method to approach this problem?

//Select Region
var roi = /* color: #98ff00 */ee.Geometry.Polygon(
        [[[-122.29786986020576, 37.4658267849554],
          [-122.07007521298897, 37.46092160855244],
          [-122.0692169061042, 37.539093572188996],
          [-122.29889982846748, 37.53882133810604]]]);

// Load MODIS land cover categories in 2001.
var landcover = ee.Image('MODIS/051/MCD12Q1/2001_01_01')
  // Select the IGBP classification band.
  .select('Land_Cover_Type_1');

Map.addLayer(landcover, {}, 'landcover');

// Load a Landsat 8 image
var image = ee.Image('LANDSAT/LC08/C01/T1_TOA/LC08_044034_20140318');
Map.setCenter(-122.1899, 37.5010, 13); // SF Bay

// Threshold the thermal band to find "hot" objects.
var hotspots = image.select('B10').gt(300);

// Mask "cold" pixels.
hotspots = hotspots.updateMask(hotspots);

// Uniquely label the patches 
var patchid = hotspots.connectedComponents(ee.Kernel.plus(1), 256);

//Add landcover band to patches and display
var img = patchid.addBands(landcover);
Map.addLayer(img.randomVisualizer(), {}, 'patch data');

// Grouped a count reducer: count number of land cover category pixels by patches
var count = img.select(['labels','Land_Cover_Type_1','B10']).reduceRegion({
  reducer: ee.Reducer.count().group({
    groupField: 0,
    groupName: 'patch id',
  }).group({
    groupField: 1,
    groupName: 'landcover',
  }),
  geometry: roi,
  scale: 30,
  maxPixels: 1e8
});

print(count);
  • Hi, how do you know that the patch id is not correct? How have you checked? May be you can illustrate that. – Rodrigo E. Principe Nov 15 '18 at 13:39
0
+200

If I understand your question correctly, then I think the following script does what you want, i.e. for every patch it shows you at which land-use category, the 50% threshold is exceeded (https://code.earthengine.google.com/58aaa2dda101c0c608c30317680055d7):

//Select Region
var roi = ee.Geometry.Polygon([[[-122.29786986020576, 37.4658267849554],
                                [-122.07007521298897, 37.46092160855244],
                                [-122.0692169061042, 37.539093572188996],
                                [-122.29889982846748, 37.53882133810604]]]);
Map.addLayer(roi);

// Load MODIS land cover categories in 2001.
var landcover = ee.Image('MODIS/051/MCD12Q1/2001_01_01').select('Land_Cover_Type_1');
Map.addLayer(landcover.clip(roi), {}, 'landcover');

// Load a Landsat 8 image
var image = ee.Image('LANDSAT/LC08/C01/T1_TOA/LC08_044034_20140318');
Map.centerObject(roi);

// Threshold the thermal band to find "hot" objects.
var hotspots = image.select('B10').gt(300);
hotspots = hotspots.updateMask(hotspots);

// Uniquely label the patches 
var patchid = hotspots.connectedComponents(ee.Kernel.plus(1), 256).select('labels');
Map.addLayer(patchid.clip(roi), {}, 'patch data');

//Add landcover band to patches and display
var img = patchid.addBands(landcover);
print(img.bandNames());

var p50 = img.reduceToVectors({
                                reducer: ee.Reducer.percentile([50]),
                                geometry: roi,
                                maxPixels: 1e15,
                              });

print(p50);
Map.addLayer(p50, {}, 'p50');

I'll focus on a single patch below to explain. If we change the reducer given to "reduceToVectors" into "ee.Reducer.histogram()", we'll get the histogram values for each patch. For example:

enter image description here

So this patch has 70 pixels with land-use value 7, 0 with value 8, 23 with value 9, ... and 229 pixels with value 13 (note the "bucketMin" of 7 and "bucketWidth" of 1). Which gives a total of 535 pixels in this specific patch.

If we sum pixels up to land-use value 10, we have 70 + 0 + 23 + 95 = 191 pixels, which is smaller than 50% of 535 = 267.5

If we now add the pixels with land-use value 12 as well, we'll have 191 + 118 = 309 pixels, which is greater than 267.5. So the land-use value at which the 50% threshold is exceeded is 12.

Changing back the reducer to the 50th percentile and looking at the same patch gives this:

enter image description here

the value of 12 is the answer you wanted (if I understood correctly).

NOTE: You mention that your patch id is not correct. For some reason, the label-value in the feature-collection (in the example 16127602200284) doesn't match with label values in the "patchid" image. I'm not sure why this is, but if that is a problem for you, you could add the correct value afterwards to the feature-collection using "reduceRegions" on the "patchid" image.

NOTE2: Add this to your code to zoom in on the example patch:

Map.setCenter(-122.17218354963791,37.49112001307097, 15)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.