[R-sig-Geo] Universal Block Kriging covariate definition for krige in gstat

Edzer Pebesma edzer.pebesma at uni-muenster.de
Wed Jan 13 15:31:03 CET 2016 


On 13/01/16 15:01, Antonio Manuel Moreno Ródenas wrote:
> Thanks a lot Edzer,
> 
> I'm not sure that would work.
> In that way I would transfer to the kriging function the averaged value
> of the covariate in the block. I'm not sure that would make the kriging
> behave correctly. 
> 
> All the points calculated with the prediction inside the block (and
> later averaged to give the block kriging prediction) will have as
> "drift" the average of the covariate in the block. Instead of getting
> the correct spatial variability inside the block (given by the
> covariate). At first sight it doesn't seems correct to me. Am I wrong?

I think so, when the operations (computing the drift, and block
averaging) are both linear, it does not matter in which order they are
carried out: f(g(x)) = g(f(x)).

It would be easy to verify by computing the universal point kriging
values and aggregating those. Try

> library(sp)
> demo(meuse, ask = FALSE, echo = FALSE)
> library(gstat)
>
> v = vgm(.5, "Sph", 900, .1)
> kr1 = krige(log(zinc)~dist, meuse, meuse.grid, v)
[using universal kriging]
>
> meuse.area$dist = aggregate(meuse.grid["dist"], meuse.area)[[1]]
> kr2 = krige(log(zinc)~dist, meuse, meuse.area, v)
[using universal kriging]
> kr2$kr1 = aggregate(kr1["var1.pred"], meuse.area)[[1]]
>
> kr2$var1.pred
[1] 5.687753
> kr2$kr1
[1] 5.685026
>
> kr2$var1.pred / kr2$kr1
[1] 1.00048
>

My guess is that the difference can be attributed to how the area is
discretized (see ?predict.gstat)

> 
> kind regards
> 
> Antonio Manuel Moreno Rodenas
> 
> /Marie Curie Early Stage Researcher/
> /PhD Candidate/
> 
> *T**U **Delft / Section Sanitary Engineering, office 4.64*____
> 
> *Civil Engineering and Geoscience Faculty *
> 
> T +31 15 278 14 62
> 
> 
> On 13 January 2016 at 14:43, Edzer Pebesma
> <edzer.pebesma at uni-muenster.de <mailto:edzer.pebesma at uni-muenster.de>>
> wrote:
> 
> 
> 
>     On 13/01/16 14:16, Antonio Manuel Moreno Ródenas wrote:
>     > Hello, I would like to rise a question on the use of predict {gstat},
>     >
>     > I'm trying to perform the estimation of a spatially distributed variable at
>     > the support scale of a particular area (Block kriging). I have access to an
>     > additional variable, it is known that the variable of interest is
>     > correlated to the new variable. So I would be interested on updating my
>     > estimation by the use of this new information. This could be done by the
>     > use of a kriging with external drift (KED), but with a block support
>     > (Universal Block kriging). Theoretically this is included in the gstat
>     > library as mentioned in the documentation.
>     >
>     > The issue comes when I try to perform the prediction:
>     >
>     > blockprediction <- predict(gstat(formula=Variabletopredict~additionalVariable,
>     > data=Observed, model=vgm), newdata = shapefile)
>     >
>     > The newdata argument should contain the prediction location. In a normal
>     > KED we would include a dataframe with a grid (coordinates in which to
>     > predict) and the values of the covariate (additionalVariable). As I'm
>     > trying to use a universal block kriging, I understood the newdata should be
>     > the region in which I'm interested to know the prediction, hence a polygon.
>     > How could I include in newdata the values of the covariate if its
>     > resolution is finer than my block?
>     >
>     > As far as I know, what block kriging does is to predict point values inside
>     > the region (which I could specified with the argument sps.args
>     > discretization), and later average them. But I don't know how to attach the
>     > covariate values to the block of interest (shapefile).
> 
>     maybe by
> 
>     shapefile = aggregate(additionalVariable, shapefile, mean)
> 
>     >
>     > Thanks in advance,
>     > I hope I could explain it properly, but I will give more details if
>     > necessary.
>     > Kind regards,
>     > Antonio
>     >
>     >       [[alternative HTML version deleted]]
>     >
>     > _______________________________________________
>     > R-sig-Geo mailing list
>     > R-sig-Geo at r-project.org <mailto: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
>     <tel:%2B49%20251%2083%2033081>
>     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 <mailto: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

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 490 bytes
Desc: OpenPGP digital signature
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20160113/9cfd9af8/attachment.bin>

More information about the R-sig-Geo mailing list