4

Does somebody know whether there is a sliding window method in R for 2d matrices and not just vectors. I need to apply median function to an image stored in matrix

Sergej Andrejev
  • 9,091
  • 11
  • 71
  • 108
  • possible duplicate of [Rolling median algorithm in C](http://stackoverflow.com/questions/1309263/rolling-median-algorithm-in-c) – Joshua Ulrich Mar 29 '12 at 19:04

2 Answers2

11

The function focal() in the excellent raster package is good for this. It takes several arguments beyond those shown in the example below, and can be used to specify a non-rectangular sliding window if that's needed.

library(raster)

## Create some example data
m <- matrix(1, ncol=10, nrow=10)
diag(m) <- 2
r <- as(m, "RasterLayer") # Coerce matrix to RasterLayer object

## Apply a function that returns a single value when passed values of cells
## in a 3-by-3 window surrounding each focal cell 
rmean <- focal(r, w=matrix(1/9, ncol=3, nrow=3), fun=mean)
rmedian <- focal(r, w=matrix(1/9, ncol=3, nrow=3), fun=median)

## Plot the results to confirm that this behaves as you'd expect
par(mfcol=c(1,3))
plot(r)
plot(rmean)
plot(rmedian)

## Coerce results back to a matrix, if you so desire
mmean <- as(rmean, "matrix")

enter image description here

Josh O'Brien
  • 159,210
  • 26
  • 366
  • 455
  • Where is focalFilter, are you sure it is in raster? – mdsumner Mar 29 '12 at 21:47
  • @mdsumner -- Oops, looks like it's now gone from the **raster** package (as explained in the **Note:** section of `?focal`). Instead, the current `focal()` can do everything that the old `focal()` and `focalFilter()` did. That's been accomplished by allowing the `w=` argument to alternatively take a matrix of weights, which is what `focalFilter()` used to do. In short, `filter()` is all that's needed any more, and I'll amend the post accordingly. Thanks for pointing that out. – Josh O'Brien Mar 29 '12 at 23:23
  • I tried it and it works. Sadly it hangs on 1280x1024 images with window width around 10. With 3 however it finishes but it takes at least 5 minutes on my computer :( – Sergej Andrejev Mar 30 '12 at 00:32
  • Too bad. I'll be interested to learn if there are any faster alternatives. I'm not too hopeful, though, as this seems like an inherently time-consuming type of calculation. (FWIW, `ma3x3.matrix()` from the **limma** package is a bit slower than `focal()`, and only allows a 3-by-3 moving window.) – Josh O'Brien Mar 30 '12 at 07:02
1

I know this is an old question, but I have come across this many times when looking to solve a similar problem. While the focal function in the raster package IS very straightforward and convienent, I have found it to be very slow when working with large rasters. There are many ways to try and address this, but one way that I found is by using system commands to "whitebox tools" which is a command-line driven set of raster analysis tools. It's main advantage it that it executes the tools in parallel and really takes advantage of multi-core CPUs. I know R has many cluster functions and packages (which I use for randomforest model raster prediction), but I have struggled with much of the parallel computing implementation in R. Whitebox tools has discrete functions for mean, max, majority, median, etc... filters (not to mention loads of terrain processing tools which is great for DEM-centric analyses).

Some example code for how I implemented a modal or majority filter (3x3 window) in R in of a large classified land cover raster (nrow=3793, ncol=6789, ncell=25750677) using whitebox tools:

system('C:/WBT2/target/release/whitebox_tools --wd="D:/Temp" ^
--run=MajorityFilter -v --input="input_rast.tif" ^
--output="maj_filt_rast.tif" --filterx=3 --filtery=3', 
wait = T, timeout=0, show.output.on.console = T)

The above code took less than 3.5 seconds to execute, meanwhile the equivalent raster package "focal" function using "modal", also from the raster package, took 5 minutes to complete coded below as:

maj_filt_rast<- focal(input_rast, fun=modal, w=matrix(1,nrow=3,ncol=3))

Getting whitebox tools compiled and installed IS a bit annoying, but good instructions are provided. In my opinion it is well worth the effort as it makes raster processes that were previously prohibitively slow in R run amazingly fast and it allows me to keep the coding for everything inside R with system commands.