123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190 |
- <h2>DESCRIPTION</h2>
- <em>v.krige</em> allows to perform kriging operations in GRASS
- environment, using R software functions in background.
- <h2>NOTES</h2>
- <em>v.krige</em> is just a front-end to R. The options and parameters
- are the same offered by packages <i>automap</i> and <i>gstat</i>.
- <p>Kriging, like other interpolation methods, is fully dependent on input
- data features. Exploratory analysis of data is encouraged to find out
- outliers, trends, anisotropies, uneven distributions and consequently
- choose the kriging algorithm that will give the most acceptable
- result. Good knowledge of the dataset is more valuable than hundreds
- of parameters or powerful hardware. See Isaaks and Srivastava's book,
- exhaustive and clear even if a bit outdated.
- <h3>Dependencies</h3>
- <dl>
- <dt><b>R software >= 2.x</b></dt>
- <dd></dd>
- <dt><b>rpy2</b></dt>
- <dd>Python binding to R. Note! <tt>rpy</tt> version 1 is not supported.</dd>
- <dt><b>R packages automap, gstat, spgrass6 and rgeos. </b></dt>
- <dd>automap is optional (provides automatic variogram fit).</dd>
- </dl>
- Install the packages via R command line (or your preferred GUI):
- <div class="code"><pre>
- install.packages("rgeos", dep=T)
- install.packages("gstat", dep=T)
- install.packages("spgrass6", dep=T)
- install.packages("automap", dep=T)
- </pre></div>
- <h4>Notes for Debian GNU/Linux</h4>
- Install the dependiencies. <b>Attention! python-rpy IS NOT
- SUITABLE.</b>:
- <div class="code"><pre>
- aptitude install R python-rpy2
- </pre></div>
- To install R packages, use either R's functions listed above (as root or as user),
- either the Debian packages [5], add to repositories' list
- for 32bit or 64bit (pick up the suitable line):
- <div class="code"><pre>
- deb <a href="http://debian.cran.r-project.org/cran2deb/debian-i386">http://debian.cran.r-project.org/cran2deb/debian-i386</a> testing/
- deb <a href="http://debian.cran.r-project.org/cran2deb/debian-amd64">http://debian.cran.r-project.org/cran2deb/debian-amd64</a> testing/
- </pre></div>
- and get the packages via aptitude:
- <div class="code"><pre>
- aptitude install r-cran-gstat r-cran-spgrass6
- </pre></div>
- <h4>Notes for Windows</h4>
- Compile GRASS following this
- <a href="http://trac.osgeo.org/grass/wiki/CompileOnWindows">guide</a>.
- You could also use Linux in a virtual machine. Or install Linux in a
- separate partition of the HD. This is not as painful as it appears,
- there are lots of guides over the Internet to help you.
- <h3>Computation time issues</h3>
- Please note that kriging calculation is slown down both by high number
- of input data points and/or high region resolution, even if they both
- contribute to a better output.
- <h2>EXAMPLES</h2>
- Kriging example based on elevation map (Spearfish data set).
- <p><b>Part 1: random sampling</b> of 2000 vector points from known
- elevation map. Each point will receive the elevation value from the
- elevation raster, as if it came from a point survey.
- <div class="code"><pre>
- g.region rast=elevation.10m -p
- v.random output=rand2k_elev n=2000
- v.db.addtable map=rand2k_elev column="elevation double precision"
- v.what.rast vect=rand2k_elev rast=elevation.10m column=elevation
- </pre></div>
- <b>Part 2: remove points lacking elevation attributes</b>. Points
- sampled at the border of the elevation map didn't receive any
- value. v.krige has no preferred action to cope with no data values, so
- the user must check for them and decide what to do (remove points,
- fill with the value of the nearest point, fill with the global/local
- mean...). In the following line of code, points with no data are
- removed from the map.
- <div class="code"><pre>
- v.extract rand2k_elev output=rand2k_elev_filt where="elevation not NULL"
- </pre></div>
- Check the result of previous line ("number of NULL attributes" must be
- 0):
- <div class="code"><pre>
- v.univar rand2k_elev_filt type=point column=elevation
- </pre></div>
- <b>Part 3: reconstruct DEM through kriging</b>. Using automatic
- variogram fit is the simplest way to run v.krige from CLI (note:
- requires R's automap package). Output map name is optional, the
- modules creates it automatically appending "_kriging" the the input
- map name and also checks for overwrite. If output_var is specified,
- the variance map is also created. Automatic variogram fit is provided
- by R package automap, the variogram models tested by the fitting
- functions are: exponential, spherical, Gaussian, Matern, M.Stein's
- parametrisation. A wider range of models is available from gstat
- package and can be tested on the GUI via the variogram plotting. If
- model is specified in the CLI, also sill, nugget and range values are
- to be provided, otherwise an error is raised (see second example of
- v.krige command).
- <div class="code"><pre>
- v.krige input=rand2k_elev_filt column=elevation output=rand2k_elev_kriging \
- output_var=rand2k_elev_kriging_var
- v.krige input=rand2k_elev_filt column=elevation output=rand2k_elev_kriging \
- output_var=rand2k_elev_kriging_var model=Lin sill=2500 nugget=0 range=1000 \
- --overwrite
- </pre></div>
- Or run wxGUI, to interactively fit the variogram and explore options:
- <div class="code"><pre>
- v.krige
- </pre></div>
- <b>Calculate prediction error</b>:
- <div class="code"><pre>
- r.mapcalc "rand2k_elev_kriging_pe = sqrt(rand2k_elev_kriging_var)"
- r.univar elevation.10m
- r.univar rand2k_elev_kriging
- r.univar rand2k_elev_kriging_pe
- </pre></div>
- The results show high errors, as the kriging techniques (ordinary and
- block kriging) are unable to handle a dataset with a trend, like the
- one used in this example: elevation is higher in the southwest corner
- and lower on northeast corner. Universal kriging can give far better
- results in these cases as it can handle the trend. It is available in
- R package gstat and will be part of a future v.krige release.
- <h2>SEE ALSO</h2>
- R package <a href="http://cran.r-project.org/web/packages/gstat/index.html">gstat</a>,
- mantained by Edzer J. Pebesma and others
- <br>
- R
- package <a href="http://cran.r-project.org/web/packages/spgrass6/index.html">spgrass6</a>,
- mantained by Roger Bivand
- <br>
- The <a href="http://grass.osgeo.org/statsgrass/grass6_r_install.html">Short
- Introduction to Geostatistical and Spatial Data Analysis with GRASS 6
- and R statistical data language</a> at the GRASS website. (includes
- installation tips)
- <br><br>
- v.krige's <a href="http://grass.osgeo.org/wiki/V.krige_GSoC_2009">wiki page</a>
- <h2>REFERENCES</h2>
- Isaaks and Srivastava, 1989: "An Introduction to Applied Geostatistics"
- (ISBN 0-19-505013-4)
- <h2>AUTHOR</h2>
- Anne Ghisla, Google Summer of Code 2009
- <p><i>Last changed: $Date$</i>
|