Using `grf`

function from R package `geoR`

, I simulated 6 replicates (each with 1000 samples) of a Gaussian random field on `[0, 1] x [0, 1]`

, with zero mean, zero nugget, and an exponential spatial covariance $$gamma(h) = expleft(frac{h}{phi}right)sigma^2,$$ with $sigma^2 = 1$ and $phi = 0.5$. In the following figure, the blue solid curve gives the semi-variogram of this true model and the blue dashed line gives the sill (which is $sigma^2 = 1$). Then for each replicate, I

- compute an empirical semi-variogram with function
`variog`

(black solid dots); - compute a monte-carlo envelope for the empirical semi-variogram (gray shaded polygon);
- fit an exponential spatial covariance model (black solid curve) to the empirical semi-variogram for $hat{sigma}^2$ and $hat{phi}$ (printed in the tile).

However, I find it difficult to understand the results / behavior as illustrated by the figure.

- I simulated data from a true model and fit that true model, why does the fitted semi-variogram differ from the true one so much? Or, why does $(hat{sigma}^2, hat{phi})$ differ from $(sigma^2, phi)$ so much? On some replicate the difference is drastic.
- How can the empirical semi-variogram (black solid dots), as well as the fitted variogram (black solid curve), breach sample variance (black dashed line)? Isn’t sample variance the sill?
- Why does the black dashed line so much away from the blue dashed line? (the answer is maybe as same as that to question 1)

Is a spatial process expected to behave like this? If I examine an AR(1) process by simulating a time series and fitting it using Yule-Walker method, I can get back to the truth, i.e., true autocorrelation coefficient and variance. This is not possible for a spatial process?

**Reproducible R code**

```
set.seed(0)
phi <- 0.5
sigmasq <- 1
## unconditional simulation of Gaussian random field with exponential covariance
EXP_GRF <- grf(1000, cov.pars = c(sigmasq, phi), cov.model = "exponential",
nsim = 6)
## estimation of semi-variogram
semi_variog <- variog(EXP_GRF)
semi_variog_i <- semi_variog
semi_variog_i$v <- NULL
par(mfrow = c(2, 3), mar = c(2.5, 2.5, 1, 1))
for (i in 1:6) {
## the i-th empirocal semi-variogram
semi_variog_i$v <- semi_variog$v[, i]
## mc envelops
env <- variog.mc.env(coords = EXP_GRF$coords, data = EXP_GRF$data[, i],
obj.variog = semi_variog_i, nsim = 99, messages = FALSE)
## plotting range
MAX <- max(env$v.lower, env$v.upper, semi_variog_i$v, sigmasq)
with(semi_variog_i, plot(u, v, ylim = c(0, MAX), type = "n"))
## empirical semi-variogram
polygon(c(semi_variog_i$u, rev(semi_variog_i$u)),
c(env$v.lower, rev(env$v.upper)), col = 8, border = NA)
with(semi_variog_i, points(u, v, pch = 19))
## variance of the data
var_data <- semi_variog$var.mark[i]
abline(h = var_data, lty = 2)
## semi-variogram of the true model
lines.variomodel(EXP_GRF, col = 4, lwd = 2)
abline(h = sigmasq, lwd = 2, col = 4, lty = 2)
## fit an variogram for parameter estimation
model <- variofit(semi_variog_i, c(var_data, phi), "exponential",
fix.nugget = TRUE, nugget = 0, message = FALSE)
## fitted variogram
lines.variomodel(model)
## add estimated `sigmasq` and `phi` as title
title(sprintf("sigmasq = %.2f, phi = %.2f", model[[2]][1], model[[2]][2]))
}
```

Cross Validated Asked by Zheyuan Li on January 7, 2022

1 AnswersI didn't have a chance to look carefully into the code, but I can already guess what you are experiencing. It may have to due with the fact that you are simulating random fields on finite domains (in a computer).

Try to increase the domain to a larger box `[0,10]x[0,10]`

or try to decrease the range of the variogram to something that is at least smaller than half of your domain size, say `0.1`

. Then plot the variograms up to a lag that is at most half of your domain. Everything that is separated by a lag that is equal or larger to half of your domain is not useful in empirical variogram calculations (too few samples, and hence spurious oscillations).

Geostatistics theory is more general than AR models, and there are theorems related to the issue you are experiencing, but I have no time at the moment to write down a more detailed answer.

Answered by juliohm on January 7, 2022

0 Asked on January 21, 2021 by atilla

1 Asked on January 21, 2021 by alajeb

1 Asked on January 21, 2021 by basketballautomation

1 Asked on January 21, 2021 by funkwecker

0 Asked on January 21, 2021

0 Asked on January 20, 2021 by igor-f

1 Asked on January 19, 2021 by wetlabstudent

1 Asked on January 19, 2021 by raghavsikaria

0 Asked on January 18, 2021 by ladan-gol

0 Asked on January 18, 2021 by thomas-moore

0 Asked on January 18, 2021 by shawn-strasser

1 Asked on January 17, 2021 by matthias

0 Asked on January 17, 2021 by cat-cuddler

case control study hypothesis testing inference observational study panel data

0 Asked on January 16, 2021 by ss-varshini

0 Asked on January 16, 2021 by sgg

gradient descent machine learning mathematical statistics risk training error

0 Asked on January 16, 2021 by adam-pollack

3 Asked on January 16, 2021 by sorcererofdm

gradient descent machine learning neural networks optimization pattern recognition

5 Asked on January 15, 2021 by aristide-herve

0 Asked on January 14, 2021 by mat

experiment design fractional factorial multivariate analysis random allocation

Get help from others!

Recent Questions

- MouseLook Script “Pops” back to the last value when the script is enabled after being disabled or destroyed
- Unity app crashes when using unmodified custom Android manifest (didn’t find class “UnityPlayerActivity”)
- How do i draw a ray in unity
- How to test consistency of responses?
- How can I understand these variograms?

Recent Answers

- eric_kernfeld on How to test consistency of responses?
- Justin Markwell on Unity app crashes when using unmodified custom Android manifest (didn’t find class “UnityPlayerActivity”)
- kjetil b halvorsen on How to test consistency of responses?
- Philipp on How do i draw a ray in unity
- DMGregory on MouseLook Script “Pops” back to the last value when the script is enabled after being disabled or destroyed

© 2022 AnswerBun.com. All rights reserved.