2

I have a raster with pixel values from 0 to 20, and I'm trying to count each pixel value with R package exactextractr. I have the raster ('def'), a SpatialPolygons object ('grid10') and an empty dataframe ('results'), pre-filled with NAs, with 21 columns x 10066 rows (there are 10066 polygons for which I want to extract counts). This is what I get:

> for (i in 0:20){
+   z <- exact_extract(def, grid10, function(value, count) {
+     count[value == i]
+   })
+   results[,i+1] <- z
+ }
  |==================================================================================================================| 100%

Error in `[<-.data.frame`(`*tmp*`, , i + 1, value = list(c(NA_real_, NA_real_,  : 
  replacement element 2 has 754 rows, need 10066
In addition: Warning message:
In `[<-.data.frame`(`*tmp*`, , i + 1, value = list(c(NA_real_, NA_real_,  :
  replacement element 1 has 10284 rows to replace 10066 rows

Can anyone suggest a correction?...I can't be far away from the correct code, as the following code works to extract the proportion coverage for each pixel value:

for (i in 0:20){
  message(i)
  z <- exact_extract(def, grid10, function(value, coverage_fraction) {
    sum(coverage_fraction[value == i],na.rm = T) / sum(coverage_fraction,na.rm = T)
  })
  results[,i+1] <- z
}

The code to set up the 'results' df was:

lenval <- length(grid10@data$left)
results <- data.frame(matrix(data=NA, nrow = length(seq(1:lenval)), 
                             ncol = length(seq(1:21))))
Improve this question

1 Answer 1

Reset to default
3

Instead of looping through 0:20 you may use table function which gives you the count per each different value.

library(exactextractr)
library(raster)
# build a raster layer
rast <- raster::raster(matrix(sample(20, 100, replace = T), ncol=10), xmn=0, ymn=0, xmx=10, ymx=10)

# and a polygon layer
poly <- sf::st_as_sfc('POLYGON ((2 2, 7 6, 4 9, 2 2))')

exact_extract(rast, poly, function(value, coverage_fraction ) {table(value)})

value [,1]
   1     4
   2     2
   3     1
   4     1
   6     1
   ...
# just to compare vs the whole raster
table(rast[])
 1  2  3  4  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 
11  4  3  7  7  4  7  4  2  1  7  7  6  3  3  4  5 10  5

# exact_extract considers also pixels that fall partially within the polygon
plot(rast); plot(poly, add = T); text(rast)

enter image description here

Improve this answer
4
  • 1
    Many thanks - with my 10066 polygons and using function(value, count) I get what looks like a list of 10066 tables - here's an extract: [[4]] < table of extent 0 > [[5]] value 2 6 7 8 9 10 11 12 13 15 16 20 1 314 270 137 139 51 6 28 1 6 7 12 [[6]] value 2 1 ....there are even a few -1 pixels I wasn't expecting! But is it possible to get this into a dataframe where the columns are pixel values (0 to 20) and the rows are the 10066 polygons?
    – Michael Smith
    Commented Jul 12, 2021 at 9:57
  • 3
    You can do something like exact_extract(rast, poly, function(df) df %>% group_by(id, value) %>% summarize(tot = sum(coverage_fraction), .groups = 'drop') %>% pivot_wider(id, names_from = value, values_from = tot), summarize_df = TRUE, include_cols = 'id') (uses dplyr and tidyr)
    – dbaston
    Commented Jul 12, 2021 at 12:44
  • 1
    Thanks @dbaston....does this mean that for exactextratr, sum * coverage fraction = count of individual values?! But it seems to work all the same.
    – Michael Smith
    Commented Jul 12, 2021 at 17:32
  • If you want to disregard the coverage fractions and count touched pixels, do tot = n() instead.
    – dbaston
    Commented Jul 12, 2021 at 23:55

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.