Average NDVI inside 95% percentile

[Leprechault]

I'd like to get NDVI average data from COPERNICUS/S2_SR collection by date using rgee package. I don like to extract the NDVI mean with all values inside my ROI, but inside the interval 5th percentile and 95th NDVI values. I try to do:

    # Packages
    library(tidyverse)
    library(rgee)
    library(sf)
    ee_Initialize(drive=TRUE)

    # Function for remove cloud and shadows ------------------------------------------
    getQABits <- function(image, qa) {
      # Convert decimal (character) to decimal (little endian)
      qa <- sum(2^(which(rev(unlist(strsplit(as.character(qa), "")) == 1))-1))
      # Return a single band image of the extracted QA bits, giving the qa value.
      image$bitwiseAnd(qa)$lt(1)
    }
    s2_clean <- function(img) {
      # Select NDVI
      img_band_selected <- img$select("B[4|8]")

      # quality band
      ndvi_qa <- img$select("QA60")

      # Select pixels to mask
      quality_mask <- getQABits(ndvi_qa, "110000000000")

      # Mask pixels with value zero.
      img_band_selected$updateMask(quality_mask)

      # Compute NDVI 
img_band_selected$select("B8")$subtract(img_band_selected$select("B4"))$divide(img_band_selected$select("B8")$add(img_band_selected$select("B4")))

     # Compute NDVI just inside 95% percentile
      img_band_selected$reduce(ee$Reducer$percentile(list(95)))
    }

    # Define a Region of interest
    roi <-ee$Geometry$Point(-52.19032,-30.25413)$buffer(500)

    # Sentinel-2 MSI dataset into the Earth Engine’s public data archive ------------              
    s2 <- ee$ImageCollection("COPERNICUS/S2_SR")

    # Select S2 images ---------------------------------------------------------------
    s2_roi  <- s2$
      filterBounds(roi)$
      filter(ee$Filter$lte("CLOUDY_PIXEL_PERCENTAGE", 1))$
      filter(ee$Filter$date(as.character(as.Date("2019-12-04")), as.character(as.Date("2020-05-03"))))$
      map(s2_clean)

    s2_roi_add_area <- s2_roi$map(
      function(img) {
        img$set("area", img$clip(roi)$geometry()$area())
      }
    )

    #Extract average NDVI values 
    ee_mean_ndvi<- ee_extract(
     x = s2_roi_add_area,
     y = roi,
     scale = 10,
     fun = ee$Reducer$mean(),
     via = "drive"
    )
    ee_mean_ndvi

I didn't find any way for combining ee$Reducer$mean() and ee$Reducer$percentile(95) . I add img_band_selected$reduce(ee$Reducer$percentile(list(95))) in s2_clean (image parameters), but I don't have changes in the NDVI results with or without this modification. Please, any help with it?

[csaybar]

Hi @Leprechault would not be better to extract the entire time series and do the filter locally?. I mean using quantile: https://www.rdocumentation.org/packages/stats/versions/3.6.2/topics/quantile.

[Leprechault]

Hi @csaybar , I don't extracts the NDVI values from the image locally because I have more than 23,000 areas since 2019 and I don't have RAM memory for aggregate this size of data after aquisition. But is strange that output contains T22JCM_p95 but the values is the same:

#  ee_mean_ndvi
#   X20191204T133221_20191204T133222_T22JCM_B8
# 1                                  0.7582231
#   X20191209T133219_20191209T133942_T22JCM_B8
# 1                                  0.7861278
#   X20191224T133221_20191224T133219_T22JCM_B8
# 1                                  0.7965397
#   X20200202T133221_20200202T133216_T22JCM_B8
# 1                                  0.7954997
#   X20200308T133219_20200308T133949_T22JCM_B8
# 1                                  0.7566285
#   X20200407T133219_20200407T133948_T22JCM_B8
# 1                                  0.8079181
#   X20200502T133231_20200502T133228_T22JCM_B8
# 1                                  0.6183093

# Average with 95% quantile NDVI data set:

#  ee_mean_ndvi
#   X20191204T133221_20191204T133222_T22JCM_p95
# 1                                   0.7582231
#   X20191209T133219_20191209T133942_T22JCM_p95
# 1                                   0.7861278
#   X20191224T133221_20191224T133219_T22JCM_p95
# 1                                   0.7965397
#   X20200202T133221_20200202T133216_T22JCM_p95
# 1                                   0.7954997
#   X20200308T133219_20200308T133949_T22JCM_p95
# 1                                   0.7566285
#   X20200407T133219_20200407T133948_T22JCM_p95
# 1                                   0.8079181
#   X20200502T133231_20200502T133228_T22JCM_p95
# 1                                   0.6183093

[csaybar]

Do you consider that ee_extract only supports reducers that only return a single value? https://r-spatial.github.io/rgee/reference/ee_extract.html

[Leprechault]

Thanks, @csaybar, I see that ee_extract using a single value using ee$Reducer$mean() or ee$Reducer$percentile(), but not both.

[csaybar]

Hi @Leprechault the problem is that ee_extract assume that the output name is the same that the reducer name. For instance, ee$Reducer$mean -> {'id': xxxxx, 'mean': '99'}

This always happens in reducers that only return a single value. To make that percentile reducer works in ee_extract run as follows:

## Not run: 
library(rgee)
library(sf)

ee_Initialize(gcs = TRUE, drive = TRUE)

# Define a Image or ImageCollection: Terraclimate
terraclimate <- ee$ImageCollection("IDAHO_EPSCOR/TERRACLIMATE") %>%
  ee$ImageCollection$filterDate("2001-01-01", "2002-01-01") %>%
  ee$ImageCollection$map(
    function(x) {
      date <- ee$Date(x$get("system:time_start"))$format('YYYY_MM_dd')
      name <- ee$String$cat("Terraclimate_pp_", date)
      x$select("pr")$rename(name)
    }
  )

# Define a geometry
nc <- st_read(
  dsn = system.file("shape/nc.shp", package = "sf"),
  stringsAsFactors = FALSE,
  quiet = TRUE
)

#Extract values - getInfo
ee_nc_rain <- ee_extract(
  x = terraclimate,
  y = nc["NAME"],
  scale = 250,
  fun = ee$Reducer$percentile(percentiles = list(95L), outputNames = list("percentile")),
  sf = TRUE
)

We are going to work in a solution to make that ee_extract works for all the EE reducers. Let me know if the problem persists!

[Leprechault]

Hi @Leprechault the problem is that ee_extract assume that the output name is the same that the reducer name. For instance, ee$Reducer$mean -> {'id': xxxxx, 'mean': '99'}

This always happens in reducers that only return a single value. To make that percentile reducer works in ee_extract run as follows:

## Not run: 
library(rgee)
library(sf)

ee_Initialize(gcs = TRUE, drive = TRUE)

# Define a Image or ImageCollection: Terraclimate
terraclimate <- ee$ImageCollection("IDAHO_EPSCOR/TERRACLIMATE") %>%
  ee$ImageCollection$filterDate("2001-01-01", "2002-01-01") %>%
  ee$ImageCollection$map(
    function(x) {
      date <- ee$Date(x$get("system:time_start"))$format('YYYY_MM_dd')
      name <- ee$String$cat("Terraclimate_pp_", date)
      x$select("pr")$rename(name)
    }
  )

# Define a geometry
nc <- st_read(
  dsn = system.file("shape/nc.shp", package = "sf"),
  stringsAsFactors = FALSE,
  quiet = TRUE
)

#Extract values - getInfo
ee_nc_rain <- ee_extract(
  x = terraclimate,
  y = nc["NAME"],
  scale = 250,
  fun = ee$Reducer$percentile(percentiles = list(95L), outputNames = list("percentile")),
  sf = TRUE
)

We are going to work in a solution to make that ee_extract works for all the EE reducers. Let me know if the problem persists!

@csaybar combine like in GEE is not a option in rgee ? (eg. var.reducers <- ee.Reducer.percentile(95L).combine(ee.Reducer.())