[STATSGRASS] gstat probability maps

[Edzer J. Pebesma]

Pierluigi,

it is a (not uncommon) misunderstanding that you would need simulation 
for this. Try:

 > v <- variogram(log(zinc)~1,loc= meuse)
 > m <- vgm(.59, "Sph", 874, .04)
 > sim <- krige(formula = log(zinc)~1, locations = meuse, model = m, 
newdata =
+ meuse.grid, nmax = 15)
[using ordinary kriging]
 > sim$p500 = 1-pnorm(log(500),sim$var1.pred,sqrt(sim$var1.var))
 > spplot(sim["p500"], main = "Pr(zn > 500)")

You could of course create a large number of simulations, count the 
fraction above 500 for each cell, but the answer would only approximate 
the one you get here, under the same model, fast and analytical.
--
Edzer

Pierluigi De Rosa wrote:
> 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 <- variogram(log(zinc)~1,loc= meuse) 
> m <- vgm(.59, "Sph", 874, .04)
> sim <- 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
>
>