3

I am trying to extract the integer index/position of the maximum value in an annual MODIS NDVI time series in Earth Engine but I am having a difficult time to figure out a plan to do this. It is quite easy to calculate the maximum value of the image collection using the following javascript code:

var bbox = ee.Geometry.Polygon([
  [-8.44299,42.20410],
  [-8.43200,41.65239],
  [-7.75634,41.65239],
  [-7.76733,42.20817],
  [-8.44299,42.20410]]);

var getNDVI = function(img){
  return img.normalizedDifference(['sur_refl_b02','sur_refl_b01']);
};

var MOD = ee.ImageCollection('MODIS/006/MOD09A1').
          filterDate('2015-01-01','2015-12-31').
          filterBounds(bbox);

var NDVI = MOD.map(getNDVI);

var NDVImax = NDVI.max();

For example, in R raster package this would correspond to an operation with calc() or stackApply() with which.max() as a function to retrieve the position of the maximum value in the raster stack time series. However, I would prefer to do this directly in GEE without having to download the data.

Suggestions?

2

A bit tricky, but you can use arrays to achieve this: https://code.earthengine.google.com/7f79588ccca7a31b9aeee7710aa1779e. For example, you can add a timestamp band to associate NDVI values with images and then sort array to find maximum values, keeping both NDVI and timestamp bands together.

// turn image collection into an array
var array = NDVI.toArray()

// sort array by the first band (NDVI)
var axes = { image:0, band:1 }
var bandIndex = 0 // index of NDVI band
var sort = array.arraySlice(axes.band, bandIndex, bandIndex + 1); 
var sorted = array.arraySort(sort);

// take the last image only (slice)
var length = sorted.arrayLength(axes.image)
var values = sorted.arraySlice(axes.image, length.subtract(1), length);

// convert back to an image  
var max = values.arrayProject([axes.band]).arrayFlatten([['ndvi', 'time']])

// visualize
Map.addLayer(max, {}, 'NDVI and time stamp', false)

var ndviMax = max.select(0)
var time = max.select(1)
Map.addLayer(ndviMax, {min: 0, max: 1}, 'NDVI' )
Map.addLayer(time.randomVisualizer(), {}, 'time stamp')

Note, that the maximum image value may vary from pixel to pixel. Also, it is better to think about image collection as a set, and not as a list, that is why I've used timestamps instead of image indices (non-existing). Image order may vary unless you sort image collection.

In case if you need to find only a single image with the overall maximum value for a given geometry, this becomes a bit easier: https://code.earthengine.google.com/2536c8d3ef38539de8274054e6c8fa3f

// compute maximum value for a given polygon for all images
NDVI = NDVI.map(function(img) {
  var max = img.reduceRegion({ reducer: ee.Reducer.max(), geometry: bbox, scale: 300 })
  return img.set({max: max.get('ndvi')})
})

// find image with a maximum value
var NDVImax = ee.Image(NDVI.sort('max', false).first())

// visualize
print('date of an image with the maximum value over bbox: ', NDVImax.date())

Map.addLayer(NDVImax, {min: 0, max: 1}, 'max')

However, in the latter case, results can be confusing, because there are multiple images with the same max value, so it might be better to use some high percentile for regional reduction instead of max.

  • Hi Gennadii! A really elegant solution for this problem. I was looking for the first case, were pixel values for the time of the maximum value vary from pixel to pixel. Definitely I have to learn more about array objects in GEE... ;-) Thanks! – Kamo Mar 26 '18 at 17:05

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.