[R-sig-Geo] apply "raster::aggregate" to multi-layer raster stack, using function that processes values from multiple layers to create a single layer?

Jonathan Greenberg jgrn at illinois.edu
Tue Jun 24 18:04:13 CEST 2014


Just to be clear, what you are really asking for is a function applied
to a focal window not at every pixel location, but centered on e.g.
every 15th pixels (so the focal windows don't overlap), correct?  If
so, as you noted this isn't what rasterEngine is designed for in the
way you are running it -- it takes 200 times longer, of course,
because each pixel using a 15x15 window is 225 pixels needed to
calculate an output.

I don't think aggregate takes more than 1 layer as an input (Robert?
Is this something that you could mod easily?) which is your problem.

You could solve this with rasterEngine is a bit of a backwards way --
if your "input" to rasterEngine (which sets the
extent/projection/resolution of the output) is actually the output,
aggregated file (blank), and one of the input variables is your input
file, you could write a function that uses the pixel position to pull
data from the input file.  projectRaster_rigorous (look at the source
code -- it is simply a rasterEngine function) shows an example of this
"start from the end" solution.  You won't want to use the polygon
cropping approach, however, you'll want to use a quicker function
(e.g. getValuesBlock would be a good one).

Basically:

input to rasterEngine: a blank aggregated file, chunk_format="raster"
function inputs: blank aggregated file, high resolution file, window size
function:
- determine the pixel coordinates of each pixel in the "chunk" of the
blank file (note that since chunk_format is set to "raster" instead of
"array", you will have access to all of the geographic information of
that chunk)
- loop through each pixel in the chunk and use its pixel coordinates
to set the parameters for getValuesBlock -- the starting row will need
to be the current pixel minus 7 pixels (if a 15x15 window)
---> special note: you will need to be careful with edge pixels, since
you may ask for a chunk that runs off the edge of your high rez image
- getValuesBlock from your high resolution file (one per pixel in the
chunk), run your analysis, return 1 values per pixel in the output
aggregated file

Does this make sense?  It is a bit counter-intuitive starting from the
end, but you should be able to get fairly rapid processing (at the
cost of increased coding) this way.  If you do follow this technique,
please follow-up with your actual completed code for the community!

--j

On Sat, Jun 21, 2014 at 11:26 AM, Carlos Carroll
<klamathconservation at gmail.com> wrote:
> While mean of calc(y, mean) would allow me to operate on the original
> stack, and then subsequently aggregate the results, what I need to do is
> operate on the entire set of values within the aggregation window (x * y *
> nlayers) in one go, because I am deriving a multivariate distance metric
> (below).
>
> multidist <- function(x,...)
> {
> footprintsize<-225
> numvars<-3
> meandist <- mean(dist(matrix(x,footprintsize,numvars), method =
> "euclidean"),na.rm=TRUE)
> return(meandist)
> }
>
>
>
> On Sat, Jun 21, 2014 at 4:35 AM, Robert J. Hijmans <r.hijmans at gmail.com>
> wrote:
>
>> Carlos,
>>
>> How about:
>>
>> x <- mean(y)
>>
>> or
>>
>> x <- calc(y, mean)
>>
>> Robert
>>
>>
>>
>> On Fri, Jun 20, 2014 at 2:50 PM, Carlos Carroll
>> <klamathconservation at gmail.com> wrote:
>> > When I use "aggregate" (package raster) on a multi-layer stack, with a
>> > function such as "mean", it returns a second multi-layer stack at the
>> > aggregated resolution. The function is applied separately to values from
>> > each layer in turn.
>> > What I want to do is send the vector/matrix of values extracted across
>> all
>> > layers to the user-defined function, which returns a single value that
>> > would then go to form the output, a single layer raster or stack.
>> > I can do this for a focal function using package spatial.tools, e.g.:
>> >
>> > v<-rasterEngine(w,fun=mean,window_dims = c(15,15,3))  ## for a 3 band
>> stack
>> > with a window size of 15 by 15.
>> >
>> > But I just need the value at the aggregate resolution, so running a focal
>> > function takes ~200 times longer than necessary in this case.
>> > Is it possible to use raster::aggregate or another function to do a
>> similar
>> > multi-layer operation?
>> >
>> >         [[alternative HTML version deleted]]
>> >
>> > _______________________________________________
>> > R-sig-Geo mailing list
>> > R-sig-Geo at r-project.org
>> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo



-- 
Jonathan A. Greenberg, PhD
Assistant Professor
Global Environmental Analysis and Remote Sensing (GEARS) Laboratory
Department of Geography and Geographic Information Science
University of Illinois at Urbana-Champaign
259 Computing Applications Building, MC-150
605 East Springfield Avenue
Champaign, IL  61820-6371
Phone: 217-300-1924
http://www.geog.illinois.edu/~jgrn/
AIM: jgrn307, MSN: jgrn307 at hotmail.com, Gchat: jgrn307, Skype: jgrn3007



More information about the R-sig-Geo mailing list