Aggregate raster by a non-integer factor with arbitrary function

Question

I would like to aggregate a population raster by a factor of 1.5, summing the values of the cells.

While aggregate() allows me to sum values when aggregating, its factor parameter accepts only integer values. projectRaster() and resample() allow me to adjust the resolution precisely, but (as far as I know) I am restricted to the prepackaged bilinear-interpolation and nearest-neighbor computation methods.

Is there a way to aggregate a raster by a non-integer factor AND specify the function to use when aggregating?

library(raster)
set.seed(10)

proj <- "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
r <- raster(resolution = 1, nrow = 100, crs = proj)
r[] <- round(rnorm(ncell(r), 100, 10))

# Doesn't accept non-integer factors
aggregate(r, fact = 1.5, fun = sum)

template <- raster(extent(r), crs = crs(r), resolution = 1.5)

# Correct resolution, but incorrect / impossible values for population
projectRaster(r, to = template, method = "ngb")
projectRaster(r, to = template, method = "bilinear")

Possible workaround

So far, the only method I've been able to come up with is to coerce the template to a SpatialPoints object; extract values from the original, higher-resolution raster; and rasterize() the result:

pts <- as(template, "SpatialPoints")
vals <- extract(r, pts)
pts2 <- SpatialPointsDataFrame(pts, data.frame(vals))

rasterize(pts2, template, field = "vals", fun = sum)

However, if the points are created at the centroids of the raster cells, I'm not sure how they are handled when extracting with a resolution of 1.5x the original raster. My preferred method would be to create a SpatialPolygonsDataFrame and rasterize with fun = mean, but (in my experience) extracting raster values using polygons is very inefficient.

Answer 1

Here a workaround:

#first resample to higher resolution
template    <- raster(extent(r), crs = crs(r), resolution = .5)  
detailedRas <- projectRaster(r, to = template, method = "ngb")

#then use an integer as a factor (in this case 3)
aggRas      <- aggregate(detailedRas, fact=3, fun=sum)

Note however that sum in this case won't return the sum of people who are living in a certain aggregated area.

I.e.: Let's say we have four cells and these values with a resolution of 1 m:

10  15  
12  18  

after resampling to 0.5 using NN:

10  10  15  15  
10  10  15  15  
12  12  18  18  
12  12  18  18  

Then aggregating by sum to 1.5m we get for the first pixel:

10+10+15+10+10+15+12+12+18 = 112

While in fact it should be something like: 10 + 15/2 + 12/2 + 18/4 = 28 (if we assume an equal population distribution over each pixel.)

I would recommend using the focal raster function with a custom / user defined function for summing up population values as you wish.

Or you divide the resampled raster by 4 and then take the sum:

2.5  2.5  3.75  3.75  
2.5  2.5  3.75  3.75   
3    3    4.5   4.5  
3    3    4.5   4.5  

2.5 + 2.5 + 3.75 + 2.5 + 2.5 + 3.75 + 3 + 3 + 4.5 = 28