具有用户定义函数的 calctest(x[1:5]) 中的 raster::calc() 错误

Question

我正在根据这个计算最大气候缺水study/scientific publication, using the Rcode that the authors publsihed with the study available here

这是一个非常简单的代码，其中输入是我从 Terraclimate 产品下载的月降雨量栅格。从每个降雨栅格（月）中减去 100 mm/month（蒸发蒸腾）后，创建一个堆栈并在 calc()

中使用以下函数

wd = stack(month.rainfall)-100 # 100 is the evapotranspiration in mm/month

# MCWD Function
mcwd.f = function(x){
result= as.numeric(x)
for(i in 1:length(result)){
wdn = result[i]
wdn1 = result[i-1]

if(i==1){
  if(wdn>0){ result[i]=0}
  else{result[i]=wdn}
}

if(i!=1){
  cwd = wdn1+wdn
  if( cwd < 0){ result[i]=cwd}
  else{result[i]=0}
   }
 }

return(result)  
 }


# Applying the Function
cwd = calc(wd, fun = mcwd.f)

但是，在我运行从 cwd 开始的行之后我得到了错误 .calcTest(x[1:5], fun, na.rm, forcefun, forceapply) 错误：不能使用这个功能

为什么会出现这个错误？我该如何解决？

我看过类似的帖子-1, 2，并使用了 na.rm=TRUE 等，但我仍然遇到同样的错误。我还看到一些帖子说 calc() 不是很好，但我不知道如何解决这个问题或找到替代方法，因为我正在使用我正在复制的研究作者发布的代码。

编辑-我使用了 100 个随机数的向量（没有 NA）来测试 mcwd.f 函数并且它有效。所以我认为这与具有 NA 像素有关，我认为如果我在 calc() 中使用 na.rm=TRUE 应该可以解决这个问题，但它并没有解决它。

Answer 1

当你问 R 问题时，总是包括一个最小的、独立的、可重现的例子。你说你的函数 "works"，但你没有提供输入和预期的输出数据。

我用你的函数得到了这个

mcwd.f(c(-5:5))
# [1]  -5  -9 -12 -14 -15 -15 -14 -12  -9  -5   0

这是你的函数的更简洁的版本

mcwd.f1 <- function(x) {
    x[1] <- min(0, x[1])
    for(i in 2:length(x)){
        x[i] <- min(0, sum(x[c(i-1, i)]))
    }
    return(x)  
}
mcwd.f1(c(-5:5))
# [1]  -5  -9 -12 -14 -15 -15 -14 -12  -9  -5   0

以及有助于解释下一个的变体

mcwd.f1a <- function(x) {
    x <- c(0, x)
    for(i in 2:length(x)){
        x[i] <- min(0, sum(x[c(i-1, i)]))
    }
    return(x[-1])  
}

矢量化版本 (after Gabor Grothendieck)

mcwd.f2 = function(x) {
    f <- function(x, y) min(x + y, 0)
    Reduce(f, x, 0, accumulate = TRUE)[-1]
}
mcwd.f2(c(-5:5))
#[1]  -5  -9 -12 -14 -15 -15 -14 -12  -9  -5   0

现在将函数与 RasterStack 和 calc

一起使用

示例数据：

library(raster)
slogo <- stack(system.file("external/rlogo.grd", package="raster")) 
s <- slogo-100

使用计算器：

x <- calc(s, mcwd.f1)
y <- calc(s, mcwd.f2)

这对这两个功能都有效（也适用于您的 mcwd.f）。尝试使用 NA 值

s[1000:3000] <- NA 
x <- calc(s, mcwd.f1)
y <- calc(s, mcwd.f2)

也有效。

x
#class      : RasterBrick 
#dimensions : 77, 101, 7777, 3  (nrow, ncol, ncell, nlayers)
#resolution : 1, 1  (x, y)
#extent     : 0, 101, 0, 77  (xmin, xmax, ymin, ymax)
#crs        : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
#source     : memory
#names      :  red, green, blue 
#min values : -100,  -200, -300 
#max values :    0,     0,    0

但不适合你的功能

y <- calc(s, mcwd.f)
Error in .calcTest(x[1:5], fun, na.rm, forcefun, forceapply) : 
  cannot use this function

当你这样做时你就会明白为什么

mcwd.f(c(-5:5,NA))
#Error in if (cwd < 0) { : missing value where TRUE/FALSE needed

为了让你的功能正常工作（但你不应该这样做，因为它效率低下）你可以像这样添加第一行

mcwd.f = function(x){
    if (any(is.na(x))) return(x*NA)
    ...

假设当存在一个 NA 时，单元格中的所有值都是或应该成为 NA

具有用户定义函数的 calctest(x[1:5]) 中的 raster::calc() 错误

raster::calc() error in calctest(x[1:5]) with userdefined function

r

raster