[GRASS-user] Searching for technical details on how v.rast.stats calculates means and other statistics

RichTC · February 20, 2015, 7:39am

Hi,

I'm trying to determine the specifics of how v.rast.stats calculates means
and other statistics, for which I have a couple of queries.

I apologise in advance if this information exists online, but I couldn't
find all the details in module guidelines or related documents and scripts.

I understand that v.rast.stats uses v.to.rast and r.univar for determining
means and statistics.

Qu. 1. Raster cells are only created (by v.to.rast) if their centre point is
within the area boundary ?
I'm using a vector area layer, and my observation after rasterizing, is that
this is the case, and cells are not created if their centre point lies
outside the area boundary.

Qu. 2. If a single cell representing a unique category in the zones layer
(generated previously by v.to.rast) overlays four cells in the input raster,
how is the mean calculated. Is it a simple sum of the four cell values
divided by n number of cells (four in this case), or is the value weighted
in some way? Not being a expert programmer, I can't really tell from the
r.univar script, and it doesn't appear from my observation to be a simple
mean calculation.

Many thanks for your comments,

Richard

GRASS version: 7.0.0svn
GRASS SVN Revision: 64042
Build Date: 2015-01-10
Build Platform: i686-pc-linux-gnu
GDAL/OGR: 1.11.1
PROJ.4: 4.9.0
GEOS: 3.4.2
SQLite: 3.7.9
Python: 2.7.3
wxPython: 2.8.12.1
Platform: Linux-3.2.0-31-generic-pae-i686-with-LinuxMint-13-maya

--
View this message in context: http://osgeo-org.1560.x6.nabble.com/Searching-for-technical-details-on-how-v-rast-stats-calculates-means-and-other-statistics-tp5188809.html
Sent from the Grass - Users mailing list archive at Nabble.com.

RichTC · February 20, 2015, 11:33am

I seem to be getting an error with r.univar when using a zone file.

I run, g.region as follows:
g.region -a raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk res=0.22

Then r.univar (see below).
There are 6 categories in the zone file, but r.univar reports null cells for
category 1? Initially I thought it was me not setting g.region correctly,
but I can't seem to get it to work as expected.

Any ideas would be welcome.

Cheers,
Richard

r.univar map=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
zones=BMR_bndry_epsg4326_rasterized@prcp05216ijk
zone 1 Category 1
total null and non-null cells: 4
total null cells: 0
Of the non-null cells:
----------------------
n: 4
minimum: 3.5821e-05
maximum: 7.74327e-05
range: 4.16117e-05
mean: 5.6472e-05
mean of absolute values: 5.6472e-05
standard deviation: 1.53325e-05
variance: 2.35085e-10
variation coefficient: 27.1506 %
sum: 0.000225888012209907
zone 2 Category 2
total null and non-null cells: 2
total null cells: 0
Of the non-null cells:
----------------------
n: 2
minimum: 3.53024e-05
maximum: 6.23036e-05
range: 2.70012e-05
mean: 4.8803e-05
mean of absolute values: 4.8803e-05
standard deviation: 1.35006e-05
variance: 1.82266e-10
variation coefficient: 27.6635 %
sum: 9.76059491222259e-05
zone 3 Category 3
total null and non-null cells: 1
total null cells: 0
Of the non-null cells:
----------------------
n: 1
minimum: 8.14342e-05
maximum: 8.14342e-05
range: 0
mean: 8.14342e-05
mean of absolute values: 8.14342e-05
standard deviation: 0
variance: 0
variation coefficient: 0 %
sum: 8.14341983641498e-05
zone 4 Category 4
total null and non-null cells: 3
total null cells: 0
Of the non-null cells:
----------------------

This is the raster file:

--
View this message in context: http://osgeo-org.1560.x6.nabble.com/Searching-for-technical-details-on-how-v-rast-stats-calculates-means-and-other-statistics-tp5188809p5188865.html
Sent from the Grass - Users mailing list archive at Nabble.com.

mlennert · February 20, 2015, 12:56pm

On 20/02/15 12:33, RichardCooper wrote:

I seem to be getting an error with r.univar when using a zone file.

I run, g.region as follows:
g.region -a raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk res=0.22

What does g.region -p show ?

Then r.univar (see below).
There are 6 categories in the zone file, but r.univar reports null cells for
category 1? Initially I thought it was me not setting g.region correctly,
but I can't seem to get it to work as expected.

What does r.univar just on each map separately show ?

Moritz

RichTC · February 20, 2015, 2:43pm

The output from g.region:
g.region -p
projection: 3 (Latitude-Longitude)
zone: 0
datum: wgs84
ellipsoid: wgs84
north: 35:38:24N
south: 15:10:48S
west: 89:45:36E
east: 140:34:48E
nsres: 0:13:12
ewres: 0:13:12
rows: 231
cols: 231
cells: 53361

This is r.univar without the zone file:
r.univar map=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
total null and non-null cells: 53361
total null cells: 16564
Of the non-null cells:
----------------------
n: 36797
minimum: 0
maximum: 0.00282167
range: 0.00282167
mean: 7.79549e-05
mean of absolute values: 7.79549e-05
standard deviation: 0.000106692
variance: 1.13832e-08
variation coefficient: 136.864 %
sum: 2.86850755307868
(Fri Feb 20 21:41:59 2015) Command finished (0 sec)

Please also see my earlier comment, where I edited the previous message:
copied here:

I'll need to check again, but there seems to be a mismatch between
categories reported in r.univar and the actual category in the zone file.
E.g. category 2 reported by r.univar is derived from category 1 of the
zones file?

--
View this message in context: http://osgeo-org.1560.x6.nabble.com/Searching-for-technical-details-on-how-v-rast-stats-calculates-means-and-other-statistics-tp5188809p5188937.html
Sent from the Grass - Users mailing list archive at Nabble.com.

Markus_Metz · February 21, 2015, 10:20am

On Fri, Feb 20, 2015 at 12:33 PM, RichardCooper <richtcooper@hotmail.com> wrote:

I seem to be getting an error with r.univar when using a zone file.

I run, g.region as follows:
g.region -a raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk res=0.22

When you use g.region -a res=, the new computational region is most
probably not identical to the raster's grid geometry. Better just use

g.region -a raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk res=0.22

If you are interested in stats for zones, you could do

g.region rast=zones
g.region align=cover

in your case

g.region rast=BMR_bndry_epsg4326_rasterized@prcp05216ijk
g.region align=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk

HTH,

Markus M

Then r.univar (see below).
There are 6 categories in the zone file, but r.univar reports null cells for
category 1? Initially I thought it was me not setting g.region correctly,
but I can't seem to get it to work as expected.

Any ideas would be welcome.

Cheers,
Richard

r.univar map=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
zones=BMR_bndry_epsg4326_rasterized@prcp05216ijk
zone 1 Category 1
total null and non-null cells: 4
total null cells: 0
Of the non-null cells:
----------------------
n: 4
minimum: 3.5821e-05
maximum: 7.74327e-05
range: 4.16117e-05
mean: 5.6472e-05
mean of absolute values: 5.6472e-05
standard deviation: 1.53325e-05
variance: 2.35085e-10
variation coefficient: 27.1506 %
sum: 0.000225888012209907
zone 2 Category 2
total null and non-null cells: 2
total null cells: 0
Of the non-null cells:
----------------------
n: 2
minimum: 3.53024e-05
maximum: 6.23036e-05
range: 2.70012e-05
mean: 4.8803e-05
mean of absolute values: 4.8803e-05
standard deviation: 1.35006e-05
variance: 1.82266e-10
variation coefficient: 27.6635 %
sum: 9.76059491222259e-05
zone 3 Category 3
total null and non-null cells: 1
total null cells: 0
Of the non-null cells:
----------------------
n: 1
minimum: 8.14342e-05
maximum: 8.14342e-05
range: 0
mean: 8.14342e-05
mean of absolute values: 8.14342e-05
standard deviation: 0
variance: 0
variation coefficient: 0 %
sum: 8.14341983641498e-05
zone 4 Category 4
total null and non-null cells: 3
total null cells: 0
Of the non-null cells:
----------------------

This is my zones file:
Type of Map: raster Number of Categories: 6
| Data Type: CELL
| Rows: 231
| Columns: 231
| Total Cells: 53361
| Projection: Latitude-Longitude
| N: 35:38:24N S: 15:10:48S Res: 0:13:12
| E: 140:34:48E W: 89:45:36E Res: 0:13:12
| Range of data: min = 1 max = 6
|
| Data Source:
| Vector Map: BMR_bndry_epsg4326@prcp05216ijk
| Original scale from vector map: 1:1
|
| Data Description:
| generated by v.to.rast
|
| Comments:
| v.to.rast input="BMR_bndry_epsg4326@prcp05216ijk" layer="1" type="ar\
| ea" output="BMR_bndry_epsg4326_rasterized" use="cat" value=1 memory=\
| 300
|

This is the raster file:

Type of Map: raster Number of Categories: 0
| Data Type: FCELL
| Rows: 230
| Columns: 230
| Total Cells: 52900
| Projection: Latitude-Longitude
| N: 35:29:24N S: 15:06:36S Res: 0:13:12
| E: 140:29:24E W: 89:53:24E Res: 0:13:12
| Range of data: min = 0 max = 0.002821672
|

--
View this message in context: http://osgeo-org.1560.x6.nabble.com/Searching-for-technical-details-on-how-v-rast-stats-calculates-means-and-other-statistics-tp5188809p5188865.html
Sent from the Grass - Users mailing list archive at Nabble.com.
_______________________________________________
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

RichTC · February 22, 2015, 6:57am

Hi Markus,

Thanks for your comments.

It still appears that v.rast.stats is giving me different values to
r.univar.

Essentially I'd like to clarify how the values are derived from
v.rast.stats. I'm guessing that the discrepancy between my r.univar and
v.rast.stats outputs relates to how I'm using g.region?

The following shows my approach for deriving statistics from both r.univar
and v.rast.stats.

Cheers,
Richard

1. Using r.univar and zones to output stats
# set g.region for rasterizing vector
g.region -p -a vector=BMR_bndry_epsg4326@prcp05216ijk res=0.22
projection: 3 (Latitude-Longitude)
zone: 0
datum: wgs84
ellipsoid: wgs84
north: 14:31:12N
south: 13:12N
west: 99:39:36E
east: 101:12E
nsres: 0:13:12
ewres: 0:13:12
rows: 6
cols: 7
cells: 42

# rasterize vector to be used for zones
v.to.rast --overwrite input=BMR_bndry_epsg4326@prcp05216ijk type=area
output=BMR_bndry_epsg4326_rasterized use=cat
Reading areas...
Reading features...
Writing raster map...
Converted areas: 6 of 0
v.to.rast complete.

# set g.region for zones with r.univar
g.region -p rast=BMR_bndry_epsg4326_rasterized@prcp05216ijk
projection: 3 (Latitude-Longitude)
zone: 0
datum: wgs84
ellipsoid: wgs84
north: 14:31:12N
south: 13:12N
west: 99:39:36E
east: 101:12E
nsres: 0:13:12
ewres: 0:13:12
rows: 6
cols: 7
cells: 42

g.region -p align=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
projection: 3 (Latitude-Longitude)
zone: 0
datum: wgs84
ellipsoid: wgs84
north: 14:35:24N
south: 13:03N
west: 99:34:12E
east: 101:19:48E
nsres: 0:13:12
ewres: 0:13:12
rows: 7
cols: 8
cells: 56

r.univar -g -t map=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
zones=BMR_bndry_epsg4326_rasterized@prcp05216ijk
zone|label|non_null_cells|null_cells|min|max|range|mean|mean_of_abs|stddev|variance|coeff_var|sum|sum_abs
1|Category
1|4|0|3.58209799742326e-05|7.74326545069925e-05|4.16116745327599e-05|5.64720030524768e-05|5.64720030524768e-05|1.53324774982512e-05|2.35084866234379e-10|27.1505820043313|0.000225888012209907|0.000225888012209907
2|Category
2|2|0|3.53023679053877e-05|6.23035812168382e-05|2.70012133114506e-05|4.88029745611129e-05|4.88029745611129e-05|1.35006066557253e-05|1.82266380072614e-10|27.6634913694027|9.76059491222259e-05|9.76059491222259e-05
3|Category
3|1|0|8.14341983641498e-05|8.14341983641498e-05|0|8.14341983641498e-05|8.14341983641498e-05|0|0|0|8.14341983641498e-05|8.14341983641498e-05
4|Category
4|3|0|5.68314899282996e-05|8.2188387750648e-05|2.53568978223484e-05|6.57776896938837e-05|6.57776896938837e-05|1.16198714234862e-05|1.35021411898352e-10|17.6653687254186|0.000197333069081651|0.000197333069081651
5|Category
5|4|0|3.15697689075023e-05|4.9253627366852e-05|1.76838584593497e-05|3.7798413359269e-05|3.7798413359269e-05|7.09251495570204e-06|5.03037683968571e-11|18.7640546926894|0.000151193653437076|0.000151193653437076
6|Category
6|1|0|1.89598704309901e-05|1.89598704309901e-05|0|1.89598704309901e-05|1.89598704309901e-05|0|0|0|1.89598704309901e-05|1.89598704309901e-05
(Sun Feb 22 13:23:02 2015) Command finished (0 sec)

2. Check outputs against v.to.rast

# set g.region for raster (stats to be derived from this layer):
g.region -p -a raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
res=0.22
projection: 3 (Latitude-Longitude)
zone: 0
datum: wgs84
ellipsoid: wgs84
north: 35:38:24N
south: 15:10:48S
west: 89:45:36E
east: 140:34:48E
nsres: 0:13:12
ewres: 0:13:12
rows: 231
cols: 231
cells: 53361

v.rast.stats map=BMR_bndry_epsg4326@prcp05216ijk
raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
column_prefix=teststat method=number,minimum,maximum,range,average,sum

Outputs from attribute table:
cat,province,teststat_number,teststat_minimum,teststat_maximum,teststat_range,teststat_average,teststat_sum
1,krung_thep_maha_nakhon_bangkok,2,3.58209799742326e-05,6.23035812168382e-05,2.64826012426056e-05,4.90622805955354e-05,9.81245611910708e-05
2,samut_prakan,1,3.53023679053877e-05,3.53023679053877e-05,0,3.53023679053877e-05,3.53023679053877e-05
3,nonthaburi,2,5.02148177474737e-05,8.14341983641498e-05,3.12193806166761e-05,6.58245080558117e-05,0.000131649016111623
4,pathum_thani,3,4.46328085672576e-05,8.2188387750648e-05,3.75555791833904e-05,6.12175620820684e-05,0.000183652686246205
5,nakhon_pathom,3,3.81464305974077e-05,6.24195599812083e-05,2.42731293838006e-05,4.99398726484893e-05,0.000149819617945468
6,samut_sakhon,2,1.89598704309901e-05,7.74326545069925e-05,5.84727840760024e-05,4.81962624689913e-05,9.63925249379827e-05

--
View this message in context: http://osgeo-org.1560.x6.nabble.com/Searching-for-technical-details-on-how-v-rast-stats-calculates-means-and-other-statistics-tp5188809p5189188.html
Sent from the Grass - Users mailing list archive at Nabble.com.

RichTC · February 22, 2015, 7:36am

I've resolved - I hadn't set g.region aligning to the raster prior to
v.to.rast:

This appears to work,with outputs from r.univar matching v.rast.stats:

g.region -p -a raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
res=0.22
g.region -p align=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
v.to.rast --overwrite input=BMR_bndry_epsg4326@prcp05216ijk type=area
output=BMR_bndry_epsg4326_rasterized use=cat
r.univar -g -t map=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
zones=BMR_bndry_epsg4326_rasterized@prcp05216ijk

v.rast.stats map=BMR_bndry_epsg4326@prcp05216ijk
raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk column_prefix=test
method=number,minimum,maximum,range,average,sum

--
View this message in context: http://osgeo-org.1560.x6.nabble.com/Searching-for-technical-details-on-how-v-rast-stats-calculates-means-and-other-statistics-tp5188809p5189189.html
Sent from the Grass - Users mailing list archive at Nabble.com.

Markus_Metz · February 24, 2015, 8:21am

On Sun, Feb 22, 2015 at 8:36 AM, RichardCooper <richtcooper@hotmail.com> wrote:

I've resolved - I hadn't set g.region aligning to the raster prior to
v.to.rast:

g.region -a res=0.22 does not align the region to the raster, but to
the resolution. This is an important difference reflected in the
results

Markus M

This appears to work,with outputs from r.univar matching v.rast.stats:

g.region -p -a raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
res=0.22
g.region -p align=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
v.to.rast --overwrite input=BMR_bndry_epsg4326@prcp05216ijk type=area
output=BMR_bndry_epsg4326_rasterized use=cat
r.univar -g -t map=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk
zones=BMR_bndry_epsg4326_rasterized@prcp05216ijk

v.rast.stats map=BMR_bndry_epsg4326@prcp05216ijk
raster=cahpa_i0apr_05216_nc_remapped_nc_1@prcp05216ijk column_prefix=test
method=number,minimum,maximum,range,average,sum

--
View this message in context: http://osgeo-org.1560.x6.nabble.com/Searching-for-technical-details-on-how-v-rast-stats-calculates-means-and-other-statistics-tp5188809p5189189.html
Sent from the Grass - Users mailing list archive at Nabble.com.
_______________________________________________
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user