<!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
<head>
  <meta content="text/html;charset=ISO-8859-1" http-equiv="Content-Type">
  <title></title>
</head>
<body bgcolor="#ffffff" text="#000000">
Dear Wolf &amp; Pierluigi, <br>
<br>
I believe that it is another common misunderstanding to believe that
indicator kriging is needed to obtain (estimates of) probabilities. You
can of course do it, try something along the lines of<br>
<pre wrap="">v &lt;- variogram(I(zinc&gt;500)~1, meuse) 
m &lt;- vgm(.15, "Sph", 874, .1) # just filling in some values here
ik &lt;- krige(I(zinc&gt;500)~1, meuse, meuse.grid, m)</pre>
but then you give up any information on how far an observation is below
or above 500, as data are transformed to 0 (below 500) and 1 (above
500) values. IMO this is valuable information, and although you may
criticize the Gaussian assumption (on the log-scale) underlying my
previous suggestion, at least it leads to "true" probabilities, instead
of interpolated 0 and 1 values that may (and will) as well lie outside
the interval [0,1].<br>
<br>
I wouldn't go as far as advising never to use indicator kriging, but
the automatic "for probabilities you need indicator kriging" clearly
misses important essentials of geostatistics.<br>
<br>
Both issues, the need to simulate or the need for indicator transforms
to obtain probabilities, have been pushed by certain authors that know
better, and re-arise from time to time. A third and related issue is
the suggestion that the (ordinary) kriging variance is not a measure
for prediction error; see the following ai-geostats thread:<br>
<a class="moz-txt-link-freetext" href="http://www.mail-archive.com/ai-geostats%40jrc.it/msg02890.html">http://www.mail-archive.com/ai-geostats%40jrc.it/msg02890.html</a><br>
<br>
Best wishes,<br>
--<br>
Edzer<br>
<br>
Wolf Bergenheim wrote:
<blockquote cite="mid46093717.9040101@bergenheim.net" type="cite">
  <pre wrap="">Pierluigi,

That would be indicator kriging:

I use gstat, so I can't give you any specific R commands, but the
principle is the same: You have to create a new dataset which is 1 if
the threshold is exceeded, and 0 if not. Now your data is 0 or 1. Then
do normal ordinary kriging on that new data column. The resulting map
shows the probability of exceeding the threshold. In gstat you have to
add  extra checks to see that your result is strictly [0.0 - 1.0]

--Wolf

On 27.03.2007 18:18, Pierluigi De Rosa wrote:
  </pre>
  <blockquote type="cite">
    <pre wrap="">Dear All,

I need to prepare a probability map of exceeding a value...
I can make the simulation, for example I tried to make it with data meuse:

v &lt;- variogram(log(zinc)~1,loc= meuse) 
m &lt;- vgm(.59, "Sph", 874, .04)
sim &lt;- krige(formula = log(zinc)~1, locations = meuse, model = m, newdata = 
meuse.grid, nmax = 15, nsim = 5)
names(sim)
[1] "sim1" "sim2" "sim3" "sim4" "sim5"

Now How can I create a probabilility map of exceding a value (for example 5)??

Thanks
Pierluigi

    </pre>
  </blockquote>
  <pre wrap=""><!---->
  </pre>
</blockquote>
<br>
</body>
</html>