Universal Block Kriging covariate definition

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
<[hidden email]<http://r-sig-geo.2731867.n2.nabble.com/user/SendEmail.jtp?type=node&node=7589380&i=0> <mailto:[hidden email]<http://r-sig-geo.2731867.n2.nabble.com/user/SendEmail.jtp?type=node&node=7589380&i=1>>>
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
    > [hidden email]<http://r-sig-geo.2731867.n2.nabble.com/user/SendEmail.jtp?type=node&node=7589380&i=2> <mailto:[hidden email]<http://r-sig-geo.2731867.n2.nabble.com/user/SendEmail.jtp?type=node&node=7589380&i=3>>
    > 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<http://www.spatialstatistics.info/>

I believe the residual variogram should then be computed using the covariate data at block support.

Sytze de Bruin
Wageningen University
Laboratory of Geo-Information Science and Remote Sensing