Cross-validation for kriging in R (package geoR): how to include the trend while reestimate is TRUE?

Previous Topic Next Topic
 
classic Classic list List threaded Threaded
3 messages Options
Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Cross-validation for kriging in R (package geoR): how to include the trend while reestimate is TRUE?

v.m.vanzoest
Dear all,

I have a question related to the function xvalid (package geoR), which I asked on StackOverflow before but unfortunately did not get answered, probably because it is too specifically related to spatial statistics and this specific function (http://stackoverflow.com/questions/43520716/cross-validation-for-kriging-in-r-how-to-include-the-trend-while-reestimating-t). I hope anyone of you is able to answer it.

I would like to compute a variogram, fit it, and then perform cross-validation. Function xvalid seems to work pretty nice to do the cross-validation. It works when I set reestimate=TRUE (so it reestimates the variogram for every point removed from the dataset in cross-validation) and it also works when using a trend. However, it does not seem to work when combining these two. Here is a reproducible example using the Meuse example dataset:

library(geoR)
library(sp)
data(meuse) # import data
coordinates(meuse) = ~x+y # make spatialpointsdataframe
meuse@proj4string <- CRS("+init=epsg:28992") # add projection
meuse_geo <- as.geodata(meuse) # create object of class geodata for geoR compatibility
meuse_geo$data <- meuse@data # attach all data (incl. covariates) to meuse_geo
meuse_vario <- variog(geodata=meuse_geo, data=meuse_geo$data$lead, trend= ~meuse_geo$data$elev) # variogram
meuse_vfit <- variofit(meuse_vario, nugget=0.1, fix.nugget=T) # fit
# cross-validation works fine:
xvalid(geodata=meuse_geo, data=meuse_geo$data$lead, model=meuse_vfit, variog.obj = meuse_vario, reestimate=F)
# cross-validation does not work when reestimate = T:
xvalid(geodata=meuse_geo, data=meuse_geo$data$lead, model=meuse_vfit, variog.obj = meuse_vario, reestimate=T)

The error I get is:

Error in variog(coords = cv.coords, data = cv.data, uvec = variog.obj$uvec,  : coords and trend have incompatible sizes

It seems to remove the point from the dataset during cross-validation, but it doesn't seem to remove the point from the covariates/trend data. Any ideas on solving this / work-arounds? Thanks a lot in advance for thinking along.

Best, Vera

---

V.M. (Vera) van Zoest, MSc | PhD candidate |
University of Twente<http://www.utwente.nl/> | Faculty ITC<http://www.itc.nl/> | Department Earth Observation Science (EOS)<https://www.itc.nl/EOS> |
ITC Building, room 2-038 | T: +31 (0)53 - 489 4412 | [hidden email]<mailto:[hidden email]> |

Study Geoinformatics: www.itc.nl/geoinformatics<http://www.itc.nl/geoinformatics>


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: Cross-validation for kriging in R (package geoR): how to include the trend while reestimate is TRUE?

Patrick Schratz
Hi Vera,

I started debugging a bit and the error is in line 192 of the `xvalid()` function which uses a subfunction `cv.f`

    res <- as.data.frame(t(apply(matrix(locations.xvalid),
      1, cv.f)))

which then does the call to `vario()` in lines 125-141.
Here, the error appears because coords (which is created in line 90) seems to be of different length then trend (which is taken from your provided variog.obj).

if (is.null(variog.obj))
            stop("xvalid: when argument reestimate = TRUE an object with the fitted variogram model must be provided in the argument variog.obj ")
          CVvar <- variog(coords = cv.coords, data = cv.data,
            uvec = variog.obj$uvec, trend = variog.obj$trend,
            lambda = variog.obj$lambda, option = variog.obj$output.type,
            estimator.type = variog.obj$estimator.type,
            nugget.tolerance = variog.obj$nugget.tolerance,
            max.dist = max(variog.obj$u), pairs.min = 2,
            bin.cloud = FALSE, direction = variog.obj$direction,
            tolerance = variog.obj$tolerance, unit.angle = "radians",
            messages = FALSE, ...)
          CVmod <- variofit(vario = CVvar, ini.cov.pars = model$cov.pars,
            cov.model = model$cov.model, fix.nugget = model$fix.nugget,
            nugget = model$nugget, fix.kappa = model$fix.kappa,
            kappa = model$kappa, max.dist = model$max.dist,
            minimisation.function = model$minimisation.function,
            weights = model$weights, messages = FALSE,
            …)

Because we are debugging somewhat deep here and the issue might be quickly solved by contacting the package authors (they should get it working quickly since they provide the option ‘reestimate = TRUE'), I would try to do so first before doing any more detailed inspection of the error.

Cheers, Patrick

PhD Student at Department of Geography - GIScience group
Friedrich-Schiller-University Jena, Germany
Tel.: +49-3641-9-48973
Web: https://pat-s.github.io

On 1. May 2017, 16:31 +0200, wrote:
>
> http://stackoverflow.com/questions/43520716/cross-validation-for-kriging-in-r-how-to-include-the-trend-while-reestimating-t

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: Cross-validation for kriging in R (package geoR): how to include the trend while reestimate is TRUE?

v.m.vanzoest
Dear Patrick,

Thanks a lot for your quick reply and taking the time for digging deeper into the function! I will contact the package authors about the issue.

Best, Vera

---

V.M. (Vera) van Zoest, MSc | PhD candidate |
University of Twente<http://www.utwente.nl/> | Faculty ITC<http://www.itc.nl/> | Department Earth Observation Science (EOS)<https://www.itc.nl/EOS> |
ITC Building, room 2-038 | T: +31 (0)53 – 489 4412 | [hidden email]<mailto:[hidden email]> |

Study Geoinformatics: www.itc.nl/geoinformatics<http://www.itc.nl/geoinformatics>


From: Patrick Schratz [mailto:[hidden email]]
Sent: maandag 1 mei 2017 20:04
To: [hidden email]; Zoest, V.M. van (ITC) <[hidden email]>
Subject: Re: [R-sig-Geo] Cross-validation for kriging in R (package geoR): how to include the trend while reestimate is TRUE?

Hi Vera,

I started debugging a bit and the error is in line 192 of the `xvalid()` function which uses a subfunction `cv.f`

    res <- as.data.frame(t(apply(matrix(locations.xvalid),
      1, cv.f)))

which then does the call to `vario()` in lines 125-141.
Here, the error appears because coords (which is created in line 90) seems to be of different length then trend (which is taken from your provided variog.obj).

if (is.null(variog.obj))
            stop("xvalid: when argument reestimate = TRUE an object with the fitted variogram model must be provided in the argument variog.obj ")
          CVvar <- variog(coords = cv.coords, data = cv.data,
            uvec = variog.obj$uvec, trend = variog.obj$trend,
            lambda = variog.obj$lambda, option = variog.obj$output.type,
            estimator.type = variog.obj$estimator.type,
            nugget.tolerance = variog.obj$nugget.tolerance,
            max.dist = max(variog.obj$u), pairs.min = 2,
            bin.cloud = FALSE, direction = variog.obj$direction,
            tolerance = variog.obj$tolerance, unit.angle = "radians",
            messages = FALSE, ...)
          CVmod <- variofit(vario = CVvar, ini.cov.pars = model$cov.pars,
            cov.model = model$cov.model, fix.nugget = model$fix.nugget,
            nugget = model$nugget, fix.kappa = model$fix.kappa,
            kappa = model$kappa, max.dist = model$max.dist,
            minimisation.function = model$minimisation.function,
            weights = model$weights, messages = FALSE,
            …)

Because we are debugging somewhat deep here and the issue might be quickly solved by contacting the package authors (they should get it working quickly since they provide the option ‘reestimate = TRUE'), I would try to do so first before doing any more detailed inspection of the error.

Cheers, Patrick

PhD Student at Department of Geography - GIScience group
Friedrich-Schiller-University Jena, Germany
Tel.: +49-3641-9-48973
Web: https://pat-s.github.io<https://pat-s.github.io/>

On 1. May 2017, 16:31 +0200, wrote:


http://stackoverflow.com/questions/43520716/cross-validation-for-kriging-in-r-how-to-include-the-trend-while-reestimating-t

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Loading...