conflict between overlay and ifelse functions in r-raster

Question

I was running a code and found a weird behaviour when using an ifelse function inside an overlay in the raster package. In short, the function will report an error if the first 5 values of every raster is NA. Why is that happening? Below is a short code that mimics the problem I found using R 3.2.3 and raster version 2.5-2, and some temporary solutions that I'm thinking on using meanwhile.

Thanks

Carlos Albero

library(raster)
cob1d <- raster(matrix(1,nr=6,nc=6))
cob1 <- cob1d; cob2 <- cob1d; cob3 <- cob1d
overlay(cob1, cob2, cob3, fun=function(x1, x2, x3)  ifelse(x1 > 0, x1 + x2 + x3,  x3))

# class       : RasterLayer 
# dimensions  : 6, 6, 36  (nrow, ncol, ncell)
# resolution  : 0.1666667, 0.1666667  (x, y)
# extent      : 0, 1, 0, 1  (xmin, xmax, ymin, ymax)
# coord. ref. : NA 
# data source : in memory
# names       : layer 
# values      : 3, 3  (min, max)


# Changing the first 5 values...

cob1[1:5] <- NA; cob2[1:5] <- NA; cob3[1:5] <- NA
overlay(cob1, cob2, cob3, fun=function(x1, x2, x3)  (x1 + x2 + x3))

gave the same result...

# but if I use `ifelse`, there is a problem:

overlay(cob1, cob2, cob3, fun=function(x1, x2, x3)  ifelse(x1 > 0, x1 + x2 + x3,  x3))
# Error in ifelse(x1, x1 + x2 + x3, x3) : 
#  argument "x2" is missing, with no default 

# Another way to solve it is adding a useless extra variable without `NA`.

cob4 <- cob1d
overlay(cob1, cob2, cob3, cob4, fun=function(x1, x2, x3, x4)  ifelse(x1 > 0, x1 + x2 + x3,  x3))
# same result as before...
# class       : RasterLayer 
# dimensions  : 6, 6, 36  (nrow, ncol, ncell)
# resolution  : 0.1666667, 0.1666667  (x, y)
# extent      : 0, 1, 0, 1  (xmin, xmax, ymin, ymax)
# coord. ref. : NA 
# data source : in memory
# names       : layer 
# values      : 3, 3  (min, max)

# or just avoiding the use of the `ifelse` function...

overlay(cob1, cob2, cob3, cob4, fun=function(x1, x2, x3, x4)  (x1 > 0) * (x1 + x2 + x3) + (x1 <= 0)*x3)

Answer 1

This is an interesting case. overlay uses the first 5 cells to determine how to process the data (via apply or via do.call). The first test is to see if apply can be used. The function used can be run via apply if all values are NA

f <- function(x1, x2, x3)  ifelse(x1 > 0, x1 + x2 + x3,  x3)
m <- matrix(NA, 5, 3)
apply(m, 1, f)

Based on that test, apply is used for all cells. However, apply fails with this function when not all values are NA:

m[1] <- 1
apply(m, 1, f)

This fails, because only has one first argument x (not three arguments as f requires).

Whereas the opposite case is common (a function that fail when only NA values are supplied) this case is rare.

raster:::calc has the argument forcefun to avoid using apply, but this is not available in overlay. I have added that to future versions of raster (version >= 2.5-4), such that you can do:

f <- function(x1, x2, x3)  ifelse(x1 > 0, x1 + x2 + x3,  x3)
overlay(cob1, cob2, cob3, fun=f, forcefun=TRUE)

In any case, you can use calc as a work-around:

s <- stack(cob1, cob2, cob3)
r <- calc(s, fun=function(x)  ifelse(x[1] > 0, x[1] + x[2] + x[3],  x[3]))

Or another alternative, based on your alternative function:

r <- (cob1 > 0) * (cob1 + cob2 + cob3) + (cob1 <= 0)*cob3