[GRASS-dev] r.series: skewness tested

Glynn Clements glynn at gclements.plus.com
Tue Jun 3 20:36:47 EDT 2008


Markus Neteler wrote:

> I have seen the r.series extension (much welcome) and tested it after
> backport with:
> 
> GRASS 6.4.svn (nc_spm_07):~ > g.region n=230000 s=214000 w=628000
> e=646000 nsres=2000 ewres=2000 -p
> 
> # numbers taken from
> http://www.nabble.com/-R--problems-with-garchFit-td7497741.html#a7501094
> i=0
> for n in `cat numbers.asc` ; do
>    i=`expr $i + 1`
>    r.mapcalc "map$i=double($n)"
> done
> r.series in=`g.mlist type=rast pat="map*" sep=','` out=skew meth=skewness
> r.info -r skew
> min=-0.639793
> max=-0.639793
> 
> R
> > mynum <- read.csv("numbers.asc", header=FALSE)
> > str(mynum)
> 'data.frame':   202 obs. of  1 variable:
>  $ V1: num   0.01071  0.02384  0.01709  0.01387 -0.00210 ...
> > skewness(mynum$V1)
> [1] -0.6371108
> 
> Is the subtle difference a precision problem somewhere?

I think that it's due to the difference between sample and population
variants of the standard deviation calculation.

I used the "population" definitions, which divide by N rather than N-1
when calculating the standard deviation. This is consistent with the
variance and standard deviation functions already provided by
lib/stats (and used by r.series).

FWIW, I get -0.641871 for /N and -0.637111 for /(N-1). The latter
matches the result from R. I'm not sure where you get -0.639793 from
(maybe there's a stray comma in the file?).

I notice that v.univar displays both the sample and population
variants of the variance and standard deviation, but only the sample
variant is given for skewness and kurtosis.

-- 
Glynn Clements <glynn at gclements.plus.com>


More information about the grass-dev mailing list