Difference between range and lag_scale

[donhalmina]

Is the lag_scale equal to the range? If so, how is it calculated because that will be wrong and must be looked into.
According to the definition of the range, it is the point(lag) beyond which the is correlation in the data (0.95*sill is reached).

I get this results when I run my data:
Stable(dim=2, var=2.63e+03, len_scale=4.05e+02, nugget=36.0, alpha=1.84)

However, when you look at the semivariogram generated the lag_scale does not correspond to what should be the range value which should be 600+ m and not 405 m. Can you kindly clarify?

[LSchueler]

There are many different definitions of the length scale for variograms.
Have a look at this example, I think that should answer all your questions.

The method you're looking for is

cov_model = Stable(dim=2, var=2.63e+03, len_scale=4.05e+02, nugget=36.0, alpha=1.84)
range = cov_model.percentile_scale(0.95)

which gives a range of around 735 m.

[donhalmina]

Thanks for the quick reply.

In that case, which length scale is used in the kriging solution and/or how can I choose which len_scale I want to use in my kriging funtion?

OrdinaryKriging(x, y, z, variogram_model=cov_model, exact_values=True)

[LSchueler]

There details are explained here.

But the easiest and also the recommended way of doing the kriging is to let GSTools estimate the variogram parameters and then pass the resulting model to the kriging.
You can find an example of how to estimate the variogram parameters here.

How did you estimate the parameters of the Stable model?

[donhalmina]

I am estimating the stable model as below.

bin_center, gamma = gs.vario_estimate((data[:, 0], data[:, 1]), data[:, 2])
nug = 36
sill = 0.95*np.var(data[:, 2])
fit_model = gs.Stable(dim=2)
fit_model.set_arg_bounds(var=[0, sill])
fit_model.set_arg_bounds(len_scale=[0, 1000])
para, pcov, r2 = fit_model.fit_variogram(bin_center, gamma, anis=False, sill = sill+nug, nugget=nug, return_r2=True)
range = fit_model.percentile_scale(0.95)
ax = fit_model.plot(x_max=max(bin_center))
ax.scatter(bin_center, gamma)
print(fit_model)
plt.show()

On the para, the len_scale is 4.05e+02 which I guess is what it using as the range and in that case, it is wrong. I say this also because, I am comparing the results obtained here with self-programmed kriging in python and the results obtained are different. The only difference/doubt with GSTools I had was how/which range value is used in estimating kriging results.

[donhalmina]

I'm reviewing your code and if I rightly understand, you're calculating the stable model as

model = nugget+ sill*(1-np.exp(-1*(lag/range)**alpha))

If i'm, right, then it should rather be (your multiple the lag by 1 instead of 3 in stable, exponential and gaussian model):

model = nugget+ sill*(1-np.exp(-3*(lag/range)**alpha))

With how you calcuate it, I get the value of range in my code = len_scale in your code.

[MuellerSeb]

If you want to have your definition of the stable variogram, you can still reimplement it as documented. There is no "right" or "wrong" with these rescaling parameters (that was the reason we introduced the rescale attribute).

import gstools as gs

class Stab(gs.CovModel):
    def default_opt_arg(self):
        return {"alpha": 1.5}

    def cor(self, h):
        return np.exp(-3 * h**self.alpha)

With Stab you get the model you wanted. Note that h will always be replaced with lag/range. This way you get the desired correlation, covariance and variogram function.

[LSchueler]

As I said, there are many different definitions of length scales for variogram models. There is no right or wrong. If you take the fitted variogram model and plug it into the krigging, than exactly that model will be used. If you prefer to work with the 95 percentile length scale, than just print it. Internally you are still using the same model.
Have a look into the link I already provided you again, you have the option for a rescaling factor. If you prefer a different kind of scaling (3 in your case), use the rescale argument in the Stable model, see here.