1

I'm trying to use gstat to predict spatially autocorrelated random variables. I'm using these sea surface temperature data. I first fit a variogram model and use this in predict. The approach works fine if I use the global dataset, but fails with some regional subsets.

My code is below:

env.dat.basin <- env.dat
# remove duplicate locations
env.dat.basin <- env.dat.basin[!duplicated(env.dat.basin[,-c(3,4)]),]

# convert to a spatial points df
coordinates(env.dat.basin) <- ~x + y
proj4string(env.dat.basin) <- CRS("+proj=longlat +datum=WGS84")

# fit variogram
v <- variogram(SST~1, env.dat.basin)
fit.v <- fit.variogram(v, model = vgm(psill = 30, model = 'Gau', range = 5000, nugget = 0), fit.sills = T, fit.ranges = T)
plot(v, fit.v)

# define a gstat object (spatial model) based on variogram fit
g.dummy <- gstat(formula=SST~1, data = env.dat.basin, beta = 1, dummy=T, model=fit.v, nmax = 30)

# make nsim simulations based on the stat object
yy <- predict(g.dummy, newdata=env.dat.basin, nsim=1)

this works, and if I subset the global dataset

env.dat.basin <- env.dat[which(env.dat$ocean == 7 | env.dat$ocean == 513), ]

it works fine too, but with a different subset

env.dat.basin <- env.dat[c(which(env.dat$ocean == 11), which(env.dat$ocean == 7 & env.dat$y <10)), ]

which looks like this:

enter image description here

the prediction returns many warnings that cause R to crash if nsim >1

Warning messages:
1: In predict.gstat(g.dummy, newdata = env.dat.basin, nsim = 1) :
  singular simulation covariance matrix
2: In predict.gstat(g.dummy, newdata = env.dat.basin, nsim = 1) :
  Covariance matrix singular at location [6.3,-21.4333,0]: skipping...
3: In predict.gstat(g.dummy, newdata = env.dat.basin, nsim = 1) :
  singular simulation covariance matrix
4: In predict.gstat(g.dummy, newdata = env.dat.basin, nsim = 1) :
  Covariance matrix singular at location [-10.6,-34.45,0]: skipping...
5: In predict.gstat(g.dummy, newdata = env.dat.basin, nsim = 1) :
  Covariance matrix singular at location [-14.885,-15.5333,0]: skipping...

I read elsewhere that this could be due to the presence of duplicate observations, but I think I removed those and the code works for the global dataset.

  • 1
    I imagine that you have no variation in the data. Have you looked at the summary statistics and how things are distributed in space (eg., using bubbleplot)? – Jeffrey Evans Nov 24 '17 at 18:59
  • @JeffreyEvans: I added a map showing the data. I'd think there's enough variation in the data, but I'd be curious to hear what you think. – Lukas Nov 25 '17 at 9:48
2

It is a consequence of using a Gaussian variogram model with zero nugget; this easily leads to near-singular covariance matrices even if points don't (exactly) overlap. Solutions: choose a different variogram model (Matern, with strong smoothing?) or add a small nugget effect (very small may already help).

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.