I want to use the spei() function in the SPEI package on a raster stack of monthly time series of water balance data. The stack has 30 years of monthly data, so 360 layers. The documentation for spei() https://cran.r-project.org/web/packages/SPEI/SPEI.pdf states the input should be "a vector, matrix or data frame with time ordered values", along with a lag (number of months). Can I run spei() without converting each pixel in the stack to a time series vector? I have explored this thread -> Apply SPEI::hargreaves function to time series from each pixel rasterbrick R which got me started setting dates as the Z dimension in the stack.
library(raster)
library(SPEI)
library(zoo)
wb.files <- list.files(path = <pathway to files>, full.names = TRUE)
wb.stack <- stack(wb.files)
dates <- seq(as.Date("1971-01-01"), as.Date("2000-12-31"), by="month")
wb.stack <- setZ(wb.stack, dates)
names(wb.stack) <- as.yearmon(getZ(wb.stack))
Then when I try to run the function, I get an error.
system.time(spei.r <- spei(wb.stack, 2))
The error I get is:
#Error in sum(is.na(data)) > 0 & na.rm == FALSE :
# operations are possible only for numeric, logical or complex types
In the thread above, the raster::overlay() function is used, and the equivalent for a single raster stack is calc() so I tried using that too.
speicalc <- function(dat, sc) {
SPEI::spei(dat, sc)
}
system.time(spei.r <- raster::calc(wb.stack, 2, fun = speicalc))
The error I get in that case is:
#Error in .calcTest(x[1:5], fun, na.rm, forcefun, forceapply) : cannot use this function
I would like to avoid collapsing every pixel to a vector or data frame, which seems like it would be very inefficient.