Skip to content
New issue

Have a question about this project? # for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “#”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? # to your account

krigeSimCE function doesn't return the simulated value when newdata has no "data" attribute #117

Open
xyzyc opened this issue Oct 17, 2022 · 2 comments

Comments

@xyzyc
Copy link

xyzyc commented Oct 17, 2022

addAttrToGeom doesn't change the newdata. I think it should be newdata = addAttrToGeom(newdata, as.data.frame(sims))
krigeSimCE = function (formula, data, newdata, model, n = 1, ext = 2)
{
stopifnot(is(model, "variogramModel"))
stopifnot(gridded(newdata))
if (!missing(data))
stopifnot(identical(data@proj4string@projargs, newdata@proj4string@projargs))
varName <- all.vars(formula[[2]])
condSim <- TRUE
if (missing(data)) {
condSim <- FALSE
message("[No data provided: performing unconditional simulation.]")
}
else {
message("[Performing conditional simulation.]")
}
covMat <- ceWrapOnTorusCalcCovRow1(newdata, model, ext = ext)
sims <- ceSim(covMat, n, newdata@grid@cells.dim, newdata@grid.index)
colnames(sims) <- paste0(varName, ".sim", 1:n)
if (!condSim) {
if ("data" %in% slotNames(newdata))
newdata@data <- cbind(newdata@data, sims)
else newdata = addAttrToGeom(newdata, as.data.frame(sims))
return(newdata)
}
obsMeanField <- krige(formula, data, newdata, model)
simMeanObsLoc <- krigeMultiple(as.formula(paste0("var1.pred ",
formula[[3]])), obsMeanField, data, model, sims)
simMeanFields <- krigeMultiple(as.formula(paste0(varName,
"
", formula[[3]])), data, newdata, model, simMeanObsLoc)
sims <- obsMeanField@data$var1.pred + sims - simMeanFields
if ("data" %in% slotNames(newdata)) {
newdata@data <- cbind(newdata@data, sims)
return(newdata)
}
newdata = addAttrToGeom(newdata, as.data.frame(sims))
}

@edzer
Copy link
Member

edzer commented Oct 17, 2022

@BenGraeler ?

@xyzyc
Copy link
Author

xyzyc commented Oct 18, 2022

I mean the "newdata = " part is missing in the original code.

# for free to join this conversation on GitHub. Already have an account? # to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants