resampled raster values out of range

Question

I would like to resample a high resolution raster to a coarser resolution, but in such a way that the maximum values of cells are retained for the coarser grid cells.

As there is no fun argument in the resample function in R's raster package, I have put together a simply custom function:

resampleCustom <- function(r1, r2) {

    resRatio <- as.integer(res(r2) / res(r1))
    ret <- aggregate(r1, fact = resRatio, fun = max)

    if (!compareRaster(ret, r2, stopiffalse = FALSE)) {
        ret <- resample(ret, r2, method = 'bilinear')
    }
    return(ret)
}

Basically, I use aggregate, where I can provide a custom function, to get close to the target raster, and then I use resample to apply some final adjustments.

I applied this to a raster that represents the projected distribution of a species of fish (where cell values represent suitability scores ranging from 0 to 1), and the odd thing is that the resulting raster has values that are greater than the max values in the original rasters.

The two rasters can be downloaded here and here.

library(raster)

# read in species raster and template
sp <- raster('Abalistes_filamentosus.tif')
template <- raster('rasterTemplate.tif')

> sp
class       : RasterLayer 
dimensions  : 360, 720, 259200  (nrow, ncol, ncell)
resolution  : 48243.14, 40790.17  (x, y)
extent      : -17367530, 17367530, -7342230, 7342230  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=cea +lon_0=0 +lat_ts=30 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs 
data source : /Users/pascaltitle/Dropbox/Abalistes_filamentosus.tif 
names       : Abalistes_filamentosus 
values      : -5.684342e-14, 1  (min, max)

> template
class       : RasterLayer 
dimensions  : 49, 116, 5684  (nrow, ncol, ncell)
resolution  : 3e+05, 3e+05  (x, y)
extent      : -17367530, 17432470, -7357770, 7342230  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=cea +lon_0=0 +lat_ts=30 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs 
data source : /Users/pascaltitle/Dropbox/rasterTemplate.tif 
names       : rasterTemplate 
values      : 1, 1  (min, max)

> resampleCustom(sp, template)
class       : RasterLayer 
dimensions  : 49, 116, 5684  (nrow, ncol, ncell)
resolution  : 3e+05, 3e+05  (x, y)
extent      : -17367530, 17432470, -7357770, 7342230  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=cea +lon_0=0 +lat_ts=30 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs 
data source : in memory
names       : Abalistes_filamentosus 
values      : -0.2061382, 1.206138  (min, max)

The max value is 1.2, but how can this be when the bilinear method should essentially be taking averages of cell values? I would expect all values of the resulting raster to be within the bounds of the original raster values.

Maybe what you are looking for is to resample using method = 'ngb' to avoid any extrapolation issues with the bilinear algorithm — dww
These are habitat suitability scores for the presence of a species, so for my purposes, if there were highly suitable cells in the fine-scale raster, I would like that to be reflected in the resulting coarse-scale raster, which is why I have been trying to implement a resampling approach that takes the max values of cells. Maybe, as you suggest, 'ngb' would be appropriate as the secondary step to resampling. — Pascal
Rather, this suggest that you should create habitat suitability scores using raster data that aligns with the low res raster that you want to aggregate it to. — Robert Hijmans

Robert Hijmans Robert Hijmans · Accepted Answer · 2018-01-22T18:14:58

The extreme values are for cells at the edge of the raster, where values are extrapolated, as there are no neighbors at one side. This shows where these values are:

x <- resampleCustom(sp, template)
a <- xyFromCell(x, which.max(x))
b <- xyFromCell(x, which.min(x))
plot(x)
points(a)
points(b)

Or

plot(Which(x < 0))
plot(Which(round(x, 15) > 0))

To remove these extreme values, you can use raster::clamp.

xc <- clamp(x, 0, 1)

By the way, what you do, first aggregate then resampling, is also what is done within raster::resample.

The fundamental problem is that your high-res raster data do not line up with the low resolution aggregation you are seeking. That suggests a mistake earlier on in your work flow. The best way to avoid this problem is probably to make the habitat suitability predictions with predictor raster data that are aligned with the high resolution raster. You perhaps did not consider that when you projected the predictor variables to +proj=cea?

resampled raster values out of range

1 Answers