v.krige.html 6.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193
  1. <h2>DESCRIPTION</h2>
  2. <em>v.krige</em> allows to perform kriging operations in GRASS
  3. environment, using R software functions in background.
  4. <h2>NOTES</h2>
  5. <em>v.krige</em> is just a front-end to R. The options and parameters
  6. are the same offered by packages <i>automap</i> and <i>gstat</i>.
  7. <p>
  8. Kriging, like other interpolation methods, is fully dependent on input
  9. data features. Exploratory analysis of data is encouraged to find out
  10. outliers, trends, anisotropies, uneven distributions and consequently
  11. choose the kriging algorithm that will give the most acceptable
  12. result. Good knowledge of the dataset is more valuable than hundreds
  13. of parameters or powerful hardware. See Isaaks and Srivastava's book,
  14. exhaustive and clear even if a bit outdated.
  15. <h3>Dependencies</h3>
  16. <dl>
  17. <dt><b>R software >= 2.x</b></dt>
  18. <dd></dd>
  19. <dt><b>rpy2</b></dt>
  20. <dd>Python binding to R. Note! <tt>rpy</tt> version 1 is not supported.</dd>
  21. <dt><b>R packages automap, gstat, spgrass6 and rgeos. </b></dt>
  22. <dd>automap is optional (provides automatic variogram fit).</dd>
  23. </dl>
  24. Install the packages via R command line (or your preferred GUI):
  25. <div class="code"><pre>
  26. install.packages("rgeos", dep=T)
  27. install.packages("gstat", dep=T)
  28. install.packages("spgrass6", dep=T)
  29. install.packages("automap", dep=T)
  30. </pre></div>
  31. <h4>Notes for Debian GNU/Linux</h4>
  32. Install the dependiencies. <b>Attention! python-rpy IS NOT
  33. SUITABLE.</b>:
  34. <div class="code"><pre>
  35. aptitude install R python-rpy2
  36. </pre></div>
  37. To install R packages, use either R's functions listed above (as root or as user),
  38. either the Debian packages [5], add to repositories' list
  39. for 32bit or 64bit (pick up the suitable line):
  40. <div class="code"><pre>
  41. deb <a href="http://debian.cran.r-project.org/cran2deb/debian-i386">http://debian.cran.r-project.org/cran2deb/debian-i386</a> testing/
  42. deb <a href="http://debian.cran.r-project.org/cran2deb/debian-amd64">http://debian.cran.r-project.org/cran2deb/debian-amd64</a> testing/
  43. </pre></div>
  44. and get the packages via aptitude:
  45. <div class="code"><pre>
  46. aptitude install r-cran-gstat r-cran-spgrass6
  47. </pre></div>
  48. <h4>Notes for Windows</h4>
  49. Compile GRASS following this
  50. <a href="http://trac.osgeo.org/grass/wiki/CompileOnWindows">guide</a>.
  51. You could also use Linux in a virtual machine. Or install Linux in a
  52. separate partition of the HD. This is not as painful as it appears,
  53. there are lots of guides over the Internet to help you.
  54. <h3>Computation time issues</h3>
  55. Please note that kriging calculation is slown down both by high number
  56. of input data points and/or high region resolution, even if they both
  57. contribute to a better output.
  58. <h2>EXAMPLES</h2>
  59. Kriging example based on elevation map (Spearfish data set).
  60. <p>
  61. <b>Part 1: random sampling</b> of 2000 vector points from known
  62. elevation map. Each point will receive the elevation value from the
  63. elevation raster, as if it came from a point survey.
  64. <div class="code"><pre>
  65. g.region rast=elevation.10m -p
  66. v.random output=rand2k_elev n=2000
  67. v.db.addtable map=rand2k_elev column="elevation double precision"
  68. v.what.rast vect=rand2k_elev rast=elevation.10m column=elevation
  69. </pre></div>
  70. <b>Part 2: remove points lacking elevation attributes</b>. Points
  71. sampled at the border of the elevation map didn't receive any
  72. value. v.krige has no preferred action to cope with no data values, so
  73. the user must check for them and decide what to do (remove points,
  74. fill with the value of the nearest point, fill with the global/local
  75. mean...). In the following line of code, points with no data are
  76. removed from the map.
  77. <div class="code"><pre>
  78. v.extract rand2k_elev output=rand2k_elev_filt where="elevation not NULL"
  79. </pre></div>
  80. Check the result of previous line ("number of NULL attributes" must be
  81. 0):
  82. <div class="code"><pre>
  83. v.univar rand2k_elev_filt type=point column=elevation
  84. </pre></div>
  85. <b>Part 3: reconstruct DEM through kriging</b>. Using automatic
  86. variogram fit is the simplest way to run v.krige from CLI (note:
  87. requires R's automap package). Output map name is optional, the
  88. modules creates it automatically appending "_kriging" the the input
  89. map name and also checks for overwrite. If output_var is specified,
  90. the variance map is also created. Automatic variogram fit is provided
  91. by R package automap, the variogram models tested by the fitting
  92. functions are: exponential, spherical, Gaussian, Matern, M.Stein's
  93. parametrisation. A wider range of models is available from gstat
  94. package and can be tested on the GUI via the variogram plotting. If
  95. model is specified in the CLI, also sill, nugget and range values are
  96. to be provided, otherwise an error is raised (see second example of
  97. v.krige command).
  98. <div class="code"><pre>
  99. v.krige input=rand2k_elev_filt column=elevation output=rand2k_elev_kriging \
  100. output_var=rand2k_elev_kriging_var
  101. v.krige input=rand2k_elev_filt column=elevation output=rand2k_elev_kriging \
  102. output_var=rand2k_elev_kriging_var model=Lin sill=2500 nugget=0 range=1000 \
  103. --overwrite
  104. </pre></div>
  105. Or run wxGUI, to interactively fit the variogram and explore options:
  106. <div class="code"><pre>
  107. v.krige
  108. </pre></div>
  109. <b>Calculate prediction error</b>:
  110. <div class="code"><pre>
  111. r.mapcalc "rand2k_elev_kriging_pe = sqrt(rand2k_elev_kriging_var)"
  112. r.univar elevation.10m
  113. r.univar rand2k_elev_kriging
  114. r.univar rand2k_elev_kriging_pe
  115. </pre></div>
  116. The results show high errors, as the kriging techniques (ordinary and
  117. block kriging) are unable to handle a dataset with a trend, like the
  118. one used in this example: elevation is higher in the southwest corner
  119. and lower on northeast corner. Universal kriging can give far better
  120. results in these cases as it can handle the trend. It is available in
  121. R package gstat and will be part of a future v.krige release.
  122. <h2>SEE ALSO</h2>
  123. R package <a href="http://cran.r-project.org/web/packages/gstat/index.html">gstat</a>,
  124. mantained by Edzer J. Pebesma and others
  125. <br>
  126. R
  127. package <a href="http://cran.r-project.org/web/packages/spgrass6/index.html">spgrass6</a>,
  128. mantained by Roger Bivand
  129. <br>
  130. The <a href="http://grass.osgeo.org/statsgrass/grass6_r_install.html">Short
  131. Introduction to Geostatistical and Spatial Data Analysis with GRASS 6
  132. and R statistical data language</a> at the GRASS website. (includes
  133. installation tips)
  134. <br><br>
  135. v.krige's <a href="http://grass.osgeo.org/wiki/V.krige_GSoC_2009">wiki page</a>
  136. <h2>REFERENCES</h2>
  137. Isaaks and Srivastava, 1989: "An Introduction to Applied Geostatistics"
  138. (ISBN 0-19-505013-4)
  139. <h2>AUTHOR</h2>
  140. Anne Ghisla, Google Summer of Code 2009
  141. <p>
  142. <i>Last changed: $Date$</i>