I am trying to calculate the Vegetation Productivity Index sensu copernicus in Google Earth Engine, which requires estimating the percentile rank. I haven't found any lead so far on how to do so, does anyone know how to do it or recommend a way to explore?

Given the example below, how would you calculate the percentile rank of the NDVI values in img compared to those in ref? What I am expecting is a value from 0-1, where the value at each pixel is the rank compared to all historical values for the same pixel.

ref = ee.ImageCollection('MODIS/006/MYD13Q1').filterDate('2006-01-01', '2006-12-31').select('NDVI') 
img = ee.ImageCollection('MODIS/006/MYD13Q1').filterDate('2007-01-01', '2007-01-10').select('NDVI').first()

2 Answers 2


By definition, Percentile Rank (not to be confused with Percentile function or Percentile) is expressed as percentage so, if you want values from 0-1 you can only divide by 100. For calculating Percentile Rank (PR), the formula is as follows; where CF—the cumulative frequency—is the count of all scores less than or equal to the score of interest, F is the frequency for the score of interest, and N is the number of scores in the distribution.

enter image description here

In your case, the score of interest is img and all scores less than or equal to the score of interest are in ref. On the other hand, MODIS/006/MYD13Q1 products must be divided by 10000 for obtaining corresponding NDVI values. I'm going to assume an arbitrary point in USA (-106.59086720561201, 38.81137001051776), however, when it is used this point with 'reduceRegion' method, some nulls values are produced with your dates range for NDVI_ref and they need to be eliminated. Following code determines PR based in above assumptions and formula, calculating its respective CF, F and N values for your date ranges.

var pt = ee.Geometry
  .Point([-106.59086720561201, 38.81137001051776]);

function divideImages(image) {
  return image.divide(10000);

var ref = ee.ImageCollection('MODIS/006/MYD13Q1')
  .filterDate('2006-01-01', '2006-12-31')

var img = ee.ImageCollection('MODIS/006/MYD13Q1')
  .filterDate('2007-01-01', '2007-01-10')

img = ee.ImageCollection(img.divide(10000));


var getNDVIvalues = function(image) {

  // Reducing region and getting value
  var ndvi_value = ee.Image(image)
    .reduceRegion(ee.Reducer.first(), pt)

  return ndvi_value;


var size1 = ref.size();

print("size list (including nulls)", size1);

var NDVI_ref = ref.toList(size1).map(getNDVIvalues);

var N = NDVI_ref.reduce(ee.Reducer.count());

print("size list (without nulls)", N);

print("NDVI ref", NDVI_ref);

var size2 = img.size();

var NDVI_img = img.toList(size2).map(getNDVIvalues).get(0);

print("NDVI img", NDVI_img);

var NDVI_ref = NDVI_ref.filter(ee.Filter.neq('item', null));

var CF_values = NDVI_ref.map(function count (ele){
  return ee.Number(ele).lte(NDVI_img).multiply(ele);

print("CF_values", CF_values);

var CF = CF_values.size();

print("CF", CF);

var F = ee.Algorithms.If(NDVI_ref.frequency(NDVI_img).eq(0), 1, NDVI_ref.frequency(NDVI_img));

print("F", F);

//PR = percentile rank
var PR = ((CF.add(ee.Number(0.5).multiply(F))).divide(N)).multiply(100);

print("PR", PR);

After running above code in GEE code editor, I got result of following image:

enter image description here

where PR of NDVI_img (0.4865) is 40.476190476190474. Complete list of NDVI_ref values are:

0.4369, 0.475, 0.3211, 0.3146, 0.204, 0.4482, 0.2881, 0.4279, 0.5854, 0.5794, null, 0.6037, 0.0748, 0.6714, 0.655, 0.6442, 0.6564, 0.7149, 0.7389, 0.709, 0.4984, 0.5799, null

They can be used for corroborating manually above formula for PR.

  • Thanks for this comprehensive answer. I've been trying to adapt this to work on full images (ee.Arrays I guess?) rather than single points but I got a bit lost. Any idea on how to achieve this?
    – e5k
    Commented Mar 31, 2021 at 2:39
  • 1
    MODIS/006/MYD13Q1 products contain global scale images (too big and too much masked pixels). So, I suggest to use your ROI for calculating the Vegetation Productivity Index placing a feature collection grid type on it. For each polygon of the grid, it would be determined its PR (mapping this feature collection; not by using an ee.Array object).
    – xunilk
    Commented Apr 1, 2021 at 2:19

I suggest doing a reduceRegion or reduceNeighborhood using ee.Reducer.percentile as your reducer.

Your code would look something like this:

var ndviPercentile = img.reduceNeighborhood({
  reducer: ee.Reducer.percentile([90]),
  geometry: inputGeometry,
  scale: /*depends on the size of your input geometry*/,
  maxPixels: 1e13
  • Percentile rank is not percentile function or percentile. You can read correct definition here: en.wikipedia.org/wiki/Percentile_rank . On the other hand, he needs to calculate PR of NDVI_img in relation to NDVI_ref values; not percentile 90.
    – xunilk
    Commented Mar 31, 2021 at 0:02
  • Ah yes, you're right, thanks!
    – M. Nicolas
    Commented Mar 31, 2021 at 18:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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