[GRASS-dev] Re: [GRASS GIS] #873: r.univar: allow percentile= to be doubles

GRASS GIS trac at osgeo.org
Wed Mar 10 12:00:47 EST 2010


#873: r.univar: allow percentile= to be doubles
--------------------------+-------------------------------------------------
  Reporter:  hamish       |       Owner:  grass-dev at lists.osgeo.org
      Type:  enhancement  |      Status:  closed                   
  Priority:  major        |   Milestone:  6.5.0                    
 Component:  Raster       |     Version:  svn-develbranch6         
Resolution:  fixed        |    Keywords:  r.univar, stats          
  Platform:  All          |         Cpu:  All                      
--------------------------+-------------------------------------------------
Changes (by hamish):

  * status:  new => closed
  * resolution:  => fixed

Comment:

 backported to relbr64 with r41377.

 Tested in 6.5svn against R's quantile(). mean residuals less than 3mm for
 Spearfish's elevation.10m.

 test method follows.

 GRASS 6.5svn:
 {{{
 #spearfish
 g.region rast=elevation.10m
 Q=`seq -s, 0 5 100`

 r.univar -e elevation.10m percentile=$Q -g
 n=2654802
 null_cells=0
 cells=2654802
 min=1061.064087
 max=1846.743408
 range=785.679321
 mean=1348.37144603811
 mean_of_abs=1348.37144603811
 stddev=175.494309916948
 variance=30798.2528132259
 coeff_var=13.0152793158443
 sum=3579659211.6848597527
 first_quartile=1196.8
 median=1309.37
 third_quartile=1480.29
 percentile_0=1061.06
 percentile_5=1123.09
 percentile_10=1153.04
 percentile_15=1171.19
 percentile_20=1184.32
 percentile_25=1196.8
 percentile_30=1210.82
 percentile_35=1229.45
 percentile_40=1251.94
 percentile_45=1280.39
 percentile_50=1309.37
 percentile_55=1346.4
 percentile_60=1378.71
 percentile_65=1410.04
 percentile_70=1440.23
 percentile_75=1480.29
 percentile_80=1525.85
 percentile_85=1568.51
 percentile_90=1613.6
 percentile_95=1671.23
 percentile_100=1846.74


 Q=`seq -s, 0 0.1 1`
 r.univar -e elevation.10m percentile=$Q -g
 n=2654802
 null_cells=0
 cells=2654802
 min=1061.064087
 max=1846.743408
 range=785.679321
 mean=1348.37144603811
 mean_of_abs=1348.37144603811
 stddev=175.494309916948
 variance=30798.2528132259
 coeff_var=13.0152793158443
 sum=3579659211.6848597527
 first_quartile=1196.8
 median=1309.37
 third_quartile=1480.29
 percentile_0=  1061.06
 percentile_0_1=1073.82
 percentile_0_2=1078.85
 percentile_0_3=1083.27
 percentile_0_4=1085.54
 percentile_0_5=1087.8
 percentile_0_6=1090.37
 percentile_0_7=1092.09
 percentile_0_8=1093.77
 percentile_0_9=1095.2
 percentile_1=  1096.48
 }}}

 R 2.7.1:
 {{{
 library(spgrass6)
 G <- gmeta6()
 spear <- rast.get6("elevation.10m", cat = FALSE, ignore.stderr = TRUE)

 summary(spear)

 Q5 <- quantile(spear$elevation.10m, probs=seq(0,1,0.05), na.rm=TRUE,
 names=TRUE, type=8)
 cat(Q5, sep="\n")
 1061.064
 1123.092
 1153.037
 1171.192
 1184.320
 1196.800
 1210.821
 1229.449
 1251.937
 1280.389
 1309.368
 1346.401
 1378.709
 1410.040
 1440.232
 1480.286
 1525.850
 1568.510
 1613.603
 1671.230
 1846.743


 Q1 <- quantile(spear$elevation.10m, probs=seq(0,0.01,0.001), na.rm=TRUE,
 names=TRUE, type=8)
 cat(Q1, sep="\n")
 1061.064
 1073.824
 1078.849
 1083.274
 1085.538
 1087.803
 1090.373
 1092.091
 1093.775
 1095.205
 1096.479
 }}}


 Octave/Matlab:
 {{{
 Q5_G = [ ...
 1061.06
 1123.09
 1153.04
 1171.19
 1184.32
 1196.8
 1210.82
 1229.45
 1251.94
 1280.39
 1309.37
 1346.4
 1378.71
 1410.04
 1440.23
 1480.29
 1525.85
 1568.51
 1613.6
 1671.23
 1846.74 ];

 Q5_R = [ ...
 1061.064
 1123.092
 1153.037
 1171.192
 1184.320
 1196.800
 1210.821
 1229.449
 1251.937
 1280.389
 1309.368
 1346.401
 1378.709
 1410.040
 1440.232
 1480.286
 1525.850
 1568.510
 1613.603
 1671.230
 1846.743 ];


 Q1_G = [ ...
 1061.06
 1073.82
 1078.85
 1083.27
 1085.54
 1087.8
 1090.37
 1092.09
 1093.77
 1095.2
 1096.48 ];

 Q1_R = [ ...
 1061.064
 1073.824
 1078.849
 1083.274
 1085.538
 1087.803
 1090.373
 1092.091
 1093.775
 1095.205
 1096.479 ];

 %summarize residuals
 >> mean ( abs(Q5_G - Q5_R) )
 ans =  0.001 meter
 >> std( abs(Q5_G - Q5_R) )
 ans = 0.0014

 >> mean ( abs(Q1_G - Q1_R) )
 ans =  0.003 meter
 >> std( abs(Q1_G - Q1_R) )
 ans = 0.0015
 }}}


 Hamish

-- 
Ticket URL: <https://trac.osgeo.org/grass/ticket/873#comment:4>
GRASS GIS <http://grass.osgeo.org>


More information about the grass-dev mailing list