Hello Bakare,

Take a look at this material from UCLA. If you're not comfortable, just tacitly ignore the math and go straight to the code section.

Cheers,

http://www.stat.ucla.edu/~nchristo/statistics_c173_c273/c173c273_lec14_w11.pdf

I see that you are mixing up between the universal kriging and the kriging with external drift. The universal kriging model is fitted when the spatial variables (easting and northing in your case) are included as covariates. So, you need a universal kriging variogram that also models the spatial variables as covariates (which you did I believe) and I would recommend you to be strict to the universal kriging model but you need to fit the CV model accordingly. You haven't set the prediction locations in the provided script, I assume it to be 'canmod.pred'. Please make sure that the prediction location also have the spatial variables in the data.

>empvarsph.var sph.mod sph.kcv yield_rk yield_rok yield_rk$var1.rk spplot(yield_rk, "var1.rk", main="RK predictions of yield")

Hope that helps.

Ordinary kriging is based on the assumption that the mean (of the random function) is an unknown constant. If the mean is not a constant then the empirical variogram does not estimate the variogram, it estimates the variogram plus the square of the mean. This is why in this case the empirical variogram will exhibit a growth rate that is quadratic or higher order. Often the mean is a function of the position coordinates and hence is estimated by fitting a trend surface (degree two or less). You need to use the residuals from the trend surface to compute the empirical variogram in that case but if you want to do kriging there are are two possibilities, one is to use the fitted variogram (fitted to the residuals) and the residual data with ordinary kriging then add back the trend surface to the kriged values. The second possibility is to use the variogram fitted to the residuals, the ORIGINAL data and universal kriging.(UK) You use the order of the polynomial of the trend surface in the UK but you do not use the coefficients of the trend surface in the UK. Do these two methods produce the same results, in general no. First of all the trend surface is not truly an estimator of the mean (note that Matheron used the term "drift" to denote the mean of the random function to distinguish it from trend surface. If you look at his 1971 lecture notes you will see that estimating the mean and fitting the variogram is a circular problem (when the mean is not constant). The UK variance incorporates the uncertainty in estimating the mean but not the uncertainty in fitting the variogram. If you add the ordinary kriging variance and the estimated error variance for the trend surface you do not get a kriging variance.

Note that it is possible that the non-constant mean is in fact a function of another variable (Matheron called this "external drift") Think of air temperature above a terrain consisting of different surface types, e.g. bare cultivated ground, paved land, land with vegetation, etc). You would first fit a regression to the original data, use the residuals to compute the empirical variogram. You have two choices just like the above. Kriging with external drift is the analogue of universal kriging.

If you want to use cross-validation then you should use UK or kriging with external drift. Interpreting the cross validation statistics is based on the fact that both UK and kriging with external drift are exact estimators (using the original data). A trend surface is not an exact predictor. i.e. the predicted value at a data location is not the same as the observed value.

Besides Universal Kriging and Kriging with External Drift you can use also Simple Kriging with Varying Loca Means (see Goovaerts, 1997, p.190), assuming that Mean is not constant but is known locally (or estimated by regression, for instance).