I've been running into all sorts of issues using ArcGIS ZonalStats and thought R could be a great way. I'm fairly new to R but I have a coding background.
The situation is that I have several rasters and a polygon shapefile with many features of different sizes (though all features are bigger than a raster cell and the polygon features are aligned to the raster). I've figured out how to get the mean value for each polygon feature using the raster library with extract:
#load packages required
require(rgdal)
require(sp)
require(raster)
# ---Set the working directory-------
datdir <- "/test_data/"
#Read in grid of water depth
ras <- raster("test_data/raster/pl_sm_rp1000/w001001.adf")
#read in polygon shape file
proxNA <- shapefile("test_data/proxy/PL_proxy_WD_NA_test")
#calc mean depth per polygon feature
#unweighted - only assigns grid to district if centroid is in that district
proxNA$RP1000 <- extract(ras, proxNA, fun = mean, na.rm = TRUE, weights = FALSE)
#plot depth values
spplot(proxNA[,'RP1000'])
The issue I have is that I also need an area based ratio between the area of the polygon and all non NA cells in the same polygon. I know what the cell size of the raster is and I can get the area for each polygon, but the missing link is the count of all non-NA cells in each feature. I managed to get the cell number of all the cells in the polygon proxNA@data$Cnumb1000 <- cellFromPolygon(ras, proxNA)
and I'm sure there is a way to get the actual value of the raster cell, which then requires a loop to get the number of all non-NA cells combined with a count, etc.
Do you have an idea or can you point me in the right direction?