[R-sig-Geo] Ordinary kriging variance of the prediction with error structure

Fri Feb 5 13:47:24 CET 2016

Thanks a lot Edzer!

Indeed, now everything goes as it should. My bad I didn't notice the
existence of the Err argument, and I was assuming that nugget would
automatically smooth the prediction.

With this everything is solved.
Best regards

On 5 February 2016 at 12:10, Edzer Pebesma <edzer.pebesma at uni-muenster.de>
wrote:

> Manuel,
>
> gstat does kriging (exact interpolatoin) when specifying a nugget
> variance, and smoothing when defining this nugget as an Err (measurement
> error) term instead:
>
> library(sp)
> library(gstat)
> demo(meuse, ask = FALSE, echo = FALSE)
> krige(log(zinc)~1, meuse, meuse[1,],vgm(.6, "Sph", 900, .05))
> krige(log(zinc)~1, meuse, meuse[1,],vgm(.6, "Sph", 900, Err = .05))
>
> # From:
>
> krige(log(zinc)~1, meuse[1:5,], meuse[1,],
>    vgm(.6, "Sph", 900, .05), set = list(debug = 32))
> krige(log(zinc)~1, meuse[1:5,], meuse[1,],
>    vgm(.6, "Sph", 900, Err = .05), set = list(debug = 32))
>
> You'll see that the difference is in the right-hand-side, the covariance
> between observations and prediction location (C0).
>
> Hth,
>
>
> On 05/02/16 11:48, Antonio Manuel Moreno Ródenas wrote:
> > Hello Ruben,
> >
> > Thanks for your answer,
> > I tried to check your reference, but still I'm having the same issue,
> >
> > Cressie 1993, he says in 3.2.27:
> > var(S|Y) = Sum(w_i)*V_io + m - c_me
> >
> > Being:
> > [w_i, vector of weights]
> > [V_io, Semivariogram value from the sampled point to the predictied]
> > [m, lagrange multiplier]
> > [c_me, nugget (due to measurement error)]
> >
> > This is equivalent to the notation I posted previously as C(h>0)= partial
> > variance + nugget - V(h>0)
> > [C, Covariance function]
> > [V, Semivariogram function]
> >
> > Substituting V_io and knowing that Sum(w_i)= 1, this gives:
> >
> > Var(S| Y) = partial variance - Sum(w_i) * C_oi - m
> > [C_oi, vector Covariance between predicted and sampled coordinates]
> >
> > This is the expression that I'm working with, in which the nugget should
> > not appear (unless x_o=x_i) isn't it?
> >
> > I'm working with the covariance matrix as the nugget will only affect
> > C(h=0) and not C(h>0) which is a property I need. And then the layout of
> > the kriging system should be the same yet adding a nugget term in the
> > diagonal of the covariance matrix.
> >
> > Am I misunderstanding something?
> >
> > Thanks and kind regards
> > Antonio
> >
> >
> > On 5 February 2016 at 10:56, rubenfcasal <rubenfcasal at gmail.com> wrote:
> >
> >> Kriging the noiseless version of Y is not “the solution of the
> >> (standard) kriging system with a nugget effect in the Covariance
> >> structure” (the semivariances/variogram at lag 0 may not be zero).
> >> See e.g. Cressie, 1993, p. 128 (for instance, eq. 3.2.27 shows the
> >> correct expression of the kriging variance).
> >>
> >> Best regards, Ruben.
> >>
> >> El 04/02/2016 a las 15:56, Antonio Manuel Moreno Ródenas escribió:
> >>> Dear r-sig-geo community,
> >>>
> >>> I would like to bring a conceptual question on the implementation of
> >>> ordinary kriging in gstat.
> >>>
> >>> I'm trying to account for my measurement error in a OK scheme. I assume
> >>> that my sampled vector Y(x_i) is a noisy realisation of S(x_i) (the
> real
> >>> variable), thus:
> >>> Y(x_i) = S(x_i) + e_i, where e_i is ~N(0,tau^2).
> >>>
> >>> If that error is assumed to follow a certain set of conditions
> (unbiased,
> >>> uncorrelated between itself/the variable and tau=constant), this is
> >>> analogous to the solution of the kriging system with a nugget effect in
> >> the
> >>> Covariance structure.
> >>>
> >>> I coded the kriging system and its solution. In order to assess if my
> >>> implementation is correct I contrasted it with the krige function in
> >> gstat.
> >>> The predicted value at each point is the same, meaning that I got
> >> correctly
> >>> the weights of the system.
> >>>
> >>> However, I'm really confused when dealing with the variance in the
> >>> prediction. Which should have this form:
> >>>
> >>> Var(S(x_o) | Y) = Var(S(x_o)) - w' * C_oi - mu
> >>>
> >>> w' weights vector
> >>> C_oi vector Covariance between predicted and sampled coordinates
> >>> mu lagrange multiplier
> >>>
> >>> If my objective would be to predict the signal of the variable S(x_o),
> >> the
> >>> term of Var(S(x_o)) will correspond to the partial variance, (*the sill
> >>> without the nugget*). This is what I understood by following the
> notation
> >>> of Model-based Geostatistics from Peter J.Diggle, where it is
> explicitly
> >>> mentioned in (pag 137 (6.8)).
> >>>
> >>> However, I only get agreement in my comparison with the krige (gstat)
> >>> variance results if I use the total sill as Var(S(x_o)) that is
> (partial
> >>> variance + nugget). So my question is:
> >>>
> >>> I'am right by thinking that still Var(S(x_o)) should not include the
> >> nugget?
> >>> What is the outcome of krige in gstat when you consider a nugget? is it
> >> the
> >>> prediction of the signal? or is it the prediction of what you would
> >> measure
> >>> at that location?
> >>>
> >>> Kind regards,
> >>>
> >>> Antonio
> >>>
> >>>       [[alternative HTML version deleted]]
> >>>
> >>> _______________________________________________
> >>> R-sig-Geo mailing list
> >>> R-sig-Geo at r-project.org
> >>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >>>
> >>
> >>
> >>         [[alternative HTML version deleted]]
> >>
> >>
> >> _______________________________________________
> >> R-sig-Geo mailing list
> >> R-sig-Geo at r-project.org
> >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >>
> >
> >       [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > R-sig-Geo at r-project.org
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >
>
> --
> Edzer Pebesma
> Institute for Geoinformatics  (ifgi),  University of Münster
> Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
> Journal of Statistical Software:   http://www.jstatsoft.org/
> Computers & Geosciences:   http://elsevier.com/locate/cageo/
> Spatial Statistics Society http://www.spatialstatistics.info
>
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

	[[alternative HTML version deleted]]