variogramfit

fits different theoretical variograms to an experimental variogram

Description

variogramfit performs a least squares fit of various theoretical variograms to an experimental, isotropic variogram. The user can choose between various bounded (e.g. spherical) and unbounded (e.g. exponential) models. A nugget variance can be modelled as well, but higher nested models are not supported.

The function works best with the function fminsearchbnd available on the FEX. You should download it from the File Exchange (File ID: #8277). If you don't have fminsearchbnd, variogramfit uses fminsearch. The problem with fminsearch is, that it might return negative variances or ranges.

Supported bounded models:
'blinear' -- bounded linear
'circular' -- circular model
'spherical' -- spherical model, =default
'pentaspherical' -- pentaspherical model

Supported unbounded functions :
'exponential' -- exponential model
'gaussian' -- gaussian variogram
'whittle' -- Whittle's elementary correlation (involves a modified Bessel function of the second kind).
'stable' -- stable models sensu Wackernagel 1995). Same as gaussian, but with different exponents. Supply the exponent alpha (<2) as an additional pn,pv-pair 'stablealpha',alpha (default = 1.5).
'matern' -- Matern model. Requires an additional pn,pv pair. 'nu',nu (shape parameter > 0, default = 1). Note that for particular values of nu the matern model reduces to other authorized variogram models. (nu = 0.5 : exponential model, nu = 1 : Whittles model, nu -> inf : Gaussian model)

References:
Wackernagel, H. (1995): Multivariate Geostatistics, Springer.
Webster, R., Oliver, M. (2001): Geostatistics for Environmental Scientists. Wiley & Sons.
Minsasny, B., McBratney, A. B. (2005): The Matérn function as general model for soil variograms. Geoderma, 3-4, 192-207.

Acknowledgements

The author wishes to acknowledge the following in the creation of this submission:
fminsearchbnd, fminsearchcon, Experimental (Semi-) Variogram
This submission has inspired the following:
Ordinary Kriging, Kriging and Inverse Distance Interpolation using GSTAT

• MATLAB release: MATLAB 7.8 (R2009a)

[Christian]

Great program Wolfgang, simplicity of use, great code commenting, everything I could have wanted...(and I didnt have to do it myself!)

[Patrick A.]

This submission is absolutely perfect in my opinion. Clear and clean code, well commented, nice code, efficient,... thanks a lot for this !

[James Ramm]

As you know I've used this in conjuction with kriging by GSTAT for large datasets.

I'm now using this to create 3D grids from scattered XYZC data ( - layered geophysical/geological data). This requires many interpolations in a loop. To this end, I wonder if it is possible to automate the model selection process so that the 'best' model may be chosen for each new sample variogram?

[Wolfgang Schwanghart]

@ James,

I think that this should be feasible without too much effort. I would, however, not want to implement it in the existing variogramfit function, but write a wrapper, that evaluates the function for various models and selects the one that has the highest coefficient of determination.

[Adam]

Please one stupid function,if I am doing right.
For this function I use output from function variogram d(val)=gammaexp and d(dist)=h?
Thank you a lot for those two functions!

[Adam]

question not function,excuse me :)

[Wolfgang Schwanghart]

Hi Adam, right!

  d.val = gammaexp
  d.dist = h

I know it's a little pedestrian...

Regards, Wolfgang

[Aditi]

Can I use this to simultaneously fit two variograms using some of the same parameter values? E.g., I want to fit both the horizontal and vertical data such that they have the same sample variance parameter.

[Wolfgang Schwanghart]

Hi Aditi,

yes, this is possible, but you have to make some changes to the function.

Note that the general variogram parameters are in vector b.

% b(1) range
% b(2) sill
% b(3) nugget

So let's say you fit your horizontal data such that
b(1) = 1;
b(2) = 5;
b(3) = []; % no nugget

Then, before fitting your next model, edit to the following lines (302-313)

% generate upper and lower bounds when fminsearchbnd is used
switch lower(params.solver)
    case {'fminsearchbnd'};
        % lower bounds
        % lb = zeros(size(b0));
        lb = [0 5];
        % upper bounds
        if nugget;
            ub = [inf max(gammaexp) max(gammaexp)]; %
        else
            % ub = [inf max(gammaexp)];
            ub = [inf 5];
        end
end

This should work. Make the obvious changes if you employ a nugget variance.

Is this a common problem and should it be implemented in variogramfit? If yes, let me know.

Best regards,
Wolfgang