Kriging returns the same value in R

Question

I have a dataframe with columns COORDS_X1, COORDS_X2, y. The dataframe has 14 observations located in a small region in Spain and I want to interpolate the unobserved locations in the region. Here is my current dataset df:

df = read.table(header = TRUE, text = "
                COORDS_X1 COORDS_X2     y
  1   280561.  4798452.  23  
2   267718.  4804594.  19  
3   270456.  4805417.  34  
4   281943.  4824197.  25  
5   283965.  4823130.  32  
6   282627.  4796964.  20  
7   286191.  4823939.  24.5
8   289349.  4794697.  23.5
9   267259.  4806210.  20  
10   274674.  4792999.  24  
11   284148.  4821649.  24  
12   273031.  4809218.  38  
13   282651.  4822399.  22  
14   281779.  4798643.  38 ")

I am no expert on this tipoc but I saw that kriging was standard methodology for this. So I am trying to use the following piece of code:

library(sp) # e.g. coordinates()
library(gstat) # e.g. variogram() 

coordinates(df) <- ~ COORDS_X1 + COORDS_X2

vgm <- variogram(y ~ 1, df)
fit <- fit.variogram(vgm, model = vgm(1, "Sph", 300, 1))

grd <- expand.grid(COORDS_X1 = seq(min(df$COORDS_X1), max(df$COORDS_X1), length.out = 100),
                   COORDS_X2 = seq(min(df$COORDS_X2), max(df$COORDS_X2), length.out = 100))
coordinates(grd) <- ~COORDS_X1 + COORDS_X2
gridded(grd) <- TRUE

krg <- krige(y ~ 1, df, grd, model = fit)

If I have a look to as.data.frame(krg), this piece of code is returning the exact same value for all the unobserved locations. Actually, if I plot the results I obtain this:

Answer 1

I suspect the issue is trying to fit a complicated model to a small number of data points in fit <- fit.variogram(vgm, model = vgm(1, "Sph", 300, 1)). Replacing this with a simpler model, e.g. fit <- fit.variogram(vgm, model = vgm("Gau")), may give better results.

From ?fit.variogram:

On singular model fits: If your variogram turns out to be a flat, horizontal or sloping line, then fitting a three parameter model such as the exponential or spherical with nugget is a bit heavy: there's an infinite number of possible combinations of sill and range (both very large) to fit to a sloping line. In this case, the returned, singular model may still be useful: just try and plot it. Gstat converges when the parameter values stabilize, and this may not be the case. Another case of singular model fits happens when a model that reaches the sill (such as the spherical) is fit with a nugget, and the range parameter starts, or converges to a value smaller than the distance of the second sample variogram estimate. In this case, again, an infinite number of possibilities occur essentially for fitting a line through a single (first sample variogram) point. In both cases, fixing one or more of the variogram model parameters may help you out.