Problems with gstat gaussian variogram and cross validation when fix nugget at 0
Stefano, Fixing the nugget at zero for the Gaussian varigoram can often cause numerical problems in kriging. If you have two observations that are close in space, their modelled semivariance will be so small that you get numerical issues with the covariance matrix. In some cases it will be singular, in your case it could be inverted, but the absolute values of the weights will be very high for some of the cross-validation locations. It seems you have two stations that are very close to each other, if you remove one of these, you can see how the results will be more similar (still not the same) to the model with a small nugget effect. Best, Jon -- Jon Olav Sk?ien European Commission Joint Research Centre ? JRC.E.1 Disaster Risk Management Unit Building 26b 1/144 | Via E.Fermi 2749, I-21027 Ispra (VA) Italy, TP 267 jon.skoien at ec.europa.eu Tel: +39 0332 789205 Disclaimer: Views expressed in this email are those of the individual and do not necessarily represent official views of the European Commission.
From: R-sig-Geo [r-sig-geo-bounces at r-project.org] on behalf of Stefano Menichetti [s.menichetti at arpat.toscana.it]
Sent: 12 October 2018 16:36
To: r-sig-geo at r-project.org
Subject: [R-sig-Geo] Problems with gstat gaussian variogram and cross validation when fix nugget at 0
Sent: 12 October 2018 16:36
To: r-sig-geo at r-project.org
Subject: [R-sig-Geo] Problems with gstat gaussian variogram and cross validation when fix nugget at 0
Dear all, I have a question about fixing nugget at 0 in a Gaussian variogram using *gstat *package I have encountered a serious problem on a variogram of piezometric data that give "strange" results if I want to fix nugget to 0. This is the file and how to import: /library(gstat)/ //file => https://drive.google.com/open?id=1AOilEhR8G8YANVuFIBUcmEg-6UG1Ajrr// /lay1 =? read.delim(paste(path,*file*,sep=""),sep="\t",dec=",")// //SpData = na.exclude(subset(lay1, select = c("X","Y","Observed"))) // //colnames(SpData)=c("X","Y","Z")/ The initial variogram: /coordinates(SpData) = ~X+Y// //var = variogram(Z~1,SpData)/ show a slow initial rise and so suggest a *Gaussian *type. After a first variogram with: /m = vgm(200,"Gau",3000,0)// // //plot(var, model=m, // //main = paste(m$model[2],"// //Nug=",m$psill[1],", Sill=",round(m$psill[2]),", Range=",round(m$range[2]),sep ="")// //)/ I have obtained a fitted variogram like this, /m = fit.variogram(var, m)/ with a reasonable cross-validation results: /krg_cv = krige.cv(Z~1, SpData, model = m)// //plot(krg_cv$observed,krg_cv$var1.pred, xlab="Z observed", ylab="Z predicted")// //abline(a=0,b=1)// //paste("media residui =",mean(krg_cv$residual))// //paste("RMS residui =",mean(krg_cv$residual^2)^0.5)/ Resulting map it is very smooth, despite the small, it seems, nugget and so *I have tried to fix nugget at 0* (exact interpolation). But cross validation results with /m = vgm(216,"Gau",3640,*0*)/ are dramatically wrong and also very different from cv results with /m = vgm(216,"Gau",3640,*1e-3*)/ What happens ? I wrong in something ? Why nugget paraneter at 0 or at small value influnce so hard the results ? Thank you very much in advance, Stefano Dott. Geol. Stefano Menichetti ============================================================= ARPAT - Agenzia per la Protezione dell'Ambiente della Toscana Settore Tecnico SIRA - Sistema Informativo Regionale Ambientale via Nicola Porpora 22 - 50144 Firenze tel. 0553206333 cell. 3668217978 (3383550147) fax 0553206410 s.menichetti at arpat.toscana.it http:\\www.arpat.toscana.it http:\\sira.arpat.toscana.it