In R, set NA cells in one raster where another raster has values

Question

I have two rasters in R. I want to set values to NA in the first raster where the second raster has values. I think this should be simple, using the raster package, with two RasterLayer objects raster1 and raster2, both the same extent and snapped to each other. They are 29775x29930.

I'm doing:

newraster <- raster1[is.na(raster2)]

But this seems to take an unnecessary amount of memory, and keeps crashing R. My computer has 8GB of memory. Is there a less memory-intensive way to do this?

Answer 1

This can easily be solved using the overlay function from the raster package. Objects rst1 and rst2 are replicates of the initial 'volcano' layer, and a random sample of n = 1000 cells in rst2 is set to NA. Afterwards, overlay is applied and the associated function rejects all cells in rst1 that hold a valid value, i.e. different from NA, in rst2.

library(raster)

rst <- raster::raster(volcano)
rst1 <- rst
rst2 <- rst

# artificial gaps
set.seed(123)
id <- sample(1:ncell(rst), 1000)
rst2[id] <- NA

# introduce na in rst1 for all locations that are non-na in rst2
rst1_mod <- overlay(rst1, rst2, fun = function(x, y) {
  x[!is.na(y[])] <- NA
  return(x)
})

plot(rst1_mod)

Answer 2

The problem is the way you are addressing the raster with newraster <- raster1[is.na(raster2)]. Try is.na(raster2) by itself to see what you get!

Trying to set values of a raster using vector terminology[] doesn't always work the way one expects. Use values() to set the actual data values of the raster:

library(raster)

# Create two small rasters and set them to random values:
r1 <- r2 <- raster(nrows=5, ncols=5)
values(r2) <- values(r1) <- rnorm(length(r1))

# Set the first row of r1 to NA:
r1[1:5] <- NA

Have a look at it:

head(r1)
           1           2          3           4           5
1         NA          NA         NA          NA          NA
2  0.1451180  0.34801689  1.0545334 -1.15284126 -0.04138286
3 -1.1370768  0.05409194 -0.7767229  0.88499661  0.34104942
4 -0.7654063 -1.03248120  1.1414939 -0.07996859 -0.34718092
5 -0.0110303 -1.70203386  0.8650742 -0.69514811 -0.31484591

Now, set the values of r2 that are not NA in r1 to NA:

values(r2)[!is.na(values(r1))] <- NA

Let's check:

head(r2)
           1         2          3        4          5
1 -0.9571857 0.4479621 -0.5601638 1.207951 -0.1459748
2         NA        NA         NA       NA         NA
3         NA        NA         NA       NA         NA
4         NA        NA         NA       NA         NA
5         NA        NA         NA       NA         NA

Answer 3

I would recommend raster::calc because according to the documentation:

For large objects calc will compute values chunk by chunk.

raster2<-reclassify(raster2,c(-Inf,Inf,0))
s->stack(raster1,raster2)
rs1<-calc(s, sum)

You might want to look into this vignette for instructions on processing large rasters.

Answer 4

This method seems like the easiest option:

s <- raster::stack(raster1, raster2)
s.new <- raster::calc(s, fun=function(x) if(sum(is.na(x)) > 0) x * NA else x)

Answer 5

Use the subs function (also possible with reclassify, calc, etc.) I demonstrate the process with a small raster first:

raster1 <- raster(ncol=99, nrow=99)
# is the raster small enough to use the following functions
canProcessInMemory(raster1)
# TRUE, so the simple method of assigning values should work
setValues(raster1, 100)
raster1[10:50, 10:50] <- NA
newraster <- subs(raster1, data.frame(old=NA, new=TRUE))

The same substitution function should work with a large raster, but you must store it on disk so that it can be processed in chunks. In my experience, any GIS will do this better than R.

raster2 <- raster(ncol=29775, nrow=29930)
canProcessInMemory(raster2)
# for such a large raster
# you will probably get the answer "FALSE" 
# meaning that setValues(raster2, 100) will not work 
# and you should write your raster to a file before working with it
# writing may take awhile
writeRaster(raster2, tempfile())