[R-sig-Geo] Calculating descriptive stats from many maps

Rainer M Krug r.m.krug at gmail.com
Tue Feb 10 09:29:25 CET 2009


Thanks a lot to all of you.

You are right Roger, I need cell-wise statistics

I like the idea of the "raster" package, and I will try it out just now.

Concerning SAGA: I'll look into that if "raster" does not work (or is to slow).

I'll report back

Rainer


On Tue, Feb 10, 2009 at 4:07 AM, Robert Hijmans <r.hijmans at gmail.com> wrote:
> Dear Rainer,
>
> This is how can you can do it with the raster package
>
> # install.packages("raster", repos="http://R-Forge.R-project.org")
> require(raster)
>
> # Try it for a few files first..
> n <- 10
>
> # create a list (or vector) of file names, e.g. :
> fn <- list()
> for (i in 1:n) { fn[i] <- paste('myfile', i, '.tif', sep='') }
>
> # make a RasterStack
> s <- stack(fn)
>
> r1 <- mCalc(s, fun=mean)
> r2 <- mCalc(s, fun=sd)
>
> #r can be plotted, coerced to sp objects, etc.
> plot(r1)
>
> # or saved to file
> r1 <- setFilename(r1, 'cellmeans.tif')
> r1 <- writeRaster(r1, format='GTiff')
>
>
> Robert
>
>
> On Tue, Feb 10, 2009 at 1:05 AM, Roger Bivand <Roger.Bivand at nhh.no> wrote:
>> On Mon, 9 Feb 2009, Tomislav Hengl wrote:
>>
>>>
>>> Dear Rainer,
>>>
>>> This is of course possible in R, and can be done in several ways:
>>>
>>> 1) for example, you can derive the average value using the rowSums
>>> function:
>>>
>>>> maps$Nsum <- rowSums(maps at data, na.rm=T, dims=1)
>>>> maps$avg <- maps$Nsum/(length(names(meuse.grid at data))-1)
>>>
>>> You could also loop the sd, mean and quantile function over a range of
>>> cells:
>>>
>>>> for(i in length(names(maps at data))) {
>>>> maps at data$sd[i] <- sd(maps at data[i,])
>>>> maps at data$mean[i] <- mean(maps at data[i,])
>>>
>>> ...
>>>>
>>>> }
>>>
>>> This could take a lot of time!
>>
>> Tom, Rainer,
>>
>> Yes, using sapply(slot(maps, "data"), summary) or similar, you get the
>> map-wise statistics. But have I misunderstood, or are the statistics in
>> question cell-wise? This would involve stacking subset areas for all 25'
>> maps, wouldn't it? Brutally, a loop in readGDAL() from rgdal with offset=
>> and region.dim= shifted? Is there a canned way to do this in the R-forge
>> raster package (by the way, regularly one of the R-forge packages showing
>> most developer activity)?
>>
>> Roger
>>
>>>
>>> 2) if your maps are rather large, try also using the SAGA function:
>>>
>>>> rsaga.get.usage(lib = "geostatistics_grid", module=5)
>>>
>>> SAGA CMD 2.0.3
>>> library path:   C:/Progra~1/saga_vc/modules
>>> library name:   geostatistics_grid
>>> module name :   Statistics for Grids
>>>
>>> This is probably the fastest method you can use.
>>>
>>> HTH
>>>
>>> T. Hengl
>>>
>>>> -----Original Message-----
>>>> From: r-sig-geo-bounces at stat.math.ethz.ch
>>>> [mailto:r-sig-geo-bounces at stat.math.ethz.ch] On Behalf
>>>> Of Rainer M Krug
>>>> Sent: Monday, February 09, 2009 4:58 PM
>>>> To: R-sig-Geo at stat.math.ethz.ch
>>>> Subject: [R-sig-Geo] Calculating descriptive stats from many maps
>>>>
>>>> Hi
>>>>
>>>> I have 25000 maps, generated by simulation predictions, covering the
>>>> same area, and would like to calculate some descriptive stats, like
>>>> mean, standard deviation, median, quartiles of all cells, to create a
>>>> "variability map".
>>>>
>>>> Is there an easy way of doing this in R?
>>>>
>>>> Thanks,
>>>>
>>>> Rainer
>>>>
>>>> --
>>>> Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation
>>>> Biology, UCT), Dipl. Phys. (Germany)
>>>>
>>>> Centre of Excellence for Invasion Biology
>>>> Faculty of Science
>>>> Natural Sciences Building
>>>> Private Bag X1
>>>> University of Stellenbosch
>>>> Matieland 7602
>>>> South Africa
>>>>
>>>> _______________________________________________
>>>> R-sig-Geo mailing list
>>>> R-sig-Geo at stat.math.ethz.ch
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo at stat.math.ethz.ch
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>>
>> --
>> Roger Bivand
>> Economic Geography Section, Department of Economics, Norwegian School of
>> Economics and Business Administration, Helleveien 30, N-5045 Bergen,
>> Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
>> e-mail: Roger.Bivand at nhh.no
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at stat.math.ethz.ch
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>



-- 
Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation
Biology, UCT), Dipl. Phys. (Germany)

Centre of Excellence for Invasion Biology
Faculty of Science
Natural Sciences Building
Private Bag X1
University of Stellenbosch
Matieland 7602
South Africa




More information about the R-sig-Geo mailing list