The best sample variogram depends on a best sampling strategy. How to optimise sampling locations see A. W. Warrick’s article (1987), “Optimization of Sampling Locations for Variogram Calculations”; Water Resources Research, 23(3): 496-500. Download paper from ResearchGate:

https://www.researchgate.net/publication/224840165_Optimization_of_sampling_locations_for_variogram_calculations

Article Optimization of sampling locations for variogram calculations

Hi Bingwei,

I'll try to answer your questions to the best of my ability.

1)   "the number Lag bins large need to than 10?  (13 in geoR ,15 in gstat by default)  Is there a standard or theory value for this?"

Re: No. There's no rule when it comes to an experimental for the same reason there's no rule when it comes to an histogram. The number of bins in an histogram or variogram (typically called lags) are user input because they're statistical tools subject to interpretation (so it's subjective in a way). Usually a good way to do it, form my point of view, is to start with a lot of lags (weaker mean - in this point there's probably to much noise in the data). Than you start to decrease the number of lags to get stronger mean (noise starts to disappear and, usually but not always, you get something that resembles a model). There are other parameters that influence a lot the experimental variogram but I'll mention them in your other questions.

2) "What is a suitable search radius of search neighborhood? How many points in the neighborhood are best case (20?)? Is this related to block kriging?"

Usually it's advised that, at the very least, your search radius includes your correlated area (the variogram ellipsoid ranges) since from that point every sample will be worth the same (sill value, the point where data is no longer correlated as a function of distance). Sometimes what is happens is that with your ranges the structures that appear in the final estimation do not make much sense (and if things get really bad you start to get image artifacts - concentric shapes in which there seems to be a failure in the smothness). You can try to correct this by increasing the number of search points or artificially increase the ranges of your variogram . A common case in 3D data (drills, logs, etc.) is to decrease the vertical search and increase the horizontal one (because vertically correlation drops quickly).

3) "The minimum sample points need large than 80 or 100? Is there a standard or theory value for this?"

It's not so much about the number of points but of the information you can get from them. With few points a variogram ellipsoid can't be conceived and thus kriging would bring no advantages. Also even if you have a lot of points but they're all aligned in one axis (for example) you won't be able to study the orthogonal. 80 to 100 would be more than enough if well distributed in space (I've done studies with much less than that). If the number of samples is a problem you could consider doing stochastic sequential simulation like sequential gaussian simulation or direct sequential simulation. These will give equiprobable images of your study. If, for instance, you would do 100 simulations and in the end you would calculate a mean image (mean for pixel in image), and variance image, you would get some that is very close to kriging (in fact it tends to the kriging image).

4) "The minimum pairs in each Lag bin need large than 3, or 10?  Is there a standard or theory value for this?"

Again is subject to interpretation. Less pairs mean less representative but this is not the only possible conclusion. Sometimes dealing with points with low number of pairs is a necessity. In fact odds say that the further the distance the more pairs exist in a point (generally speaking). The first thing I would do is ask myself if that point with almost no pairs makes sense with the general model the variogram seems to have or not. If yes no problem. If no I would have to try to make it more representative (perhaps by decreasing the number of lags and increasing the size of each, or by increasing the tolerance). Sometimes  points with few pair are local phenomena and typically you are trying to get a global model of anisotropy. Remember it's statistics, no right answer, just right interpretations.

5) "The first several points in sample variogram most important,  Is there a standard value or theory to know how many points in the first  part of lag distance?"

They're are the most important in most cases but not always. Since you're trying, I believe, to do kriging estimation te first questions you should ask yourself is: in most nodes of my grid what are the most likely distances to be used for each node? Imagine that you have samples that are more or less 100 meters apart from each other. Most nodes will be within 30 to 70 meters from any sample. Shouldn't these have special attention? It all depends on the case you're studying, but remember the final purpose of variogram modelling is to build a variogram ellipsoid that will give a satisfactory result in the estimation. So the question remains: From the point of view of the kriging algorithm is the variogram ellipsoid you've chosen a good one?

6) "What is a suitable range (the average size of physical anomalies in the

spatial fluctuation?)and sill threshold for the data sets? Can we do some prediction or estimation of them before draw the sample variogram based on knowledge or experience."

From my point of view you can and you should. If the final estimation does not make sense for the structures you're trying to conceive than probably something is wrong. Sometimes you even have to artificially increase or decrease the ranges to be able to have a more correct "version" of the estimation. It's heuristics actually. It all depends on your experience and knowledge of the phenomena you are trying to study (just as you mentioned). Sometimes when my students do kriging estimations they get these sphere like structures all around the map (it's usually because of the variograms) and I ask them: Have you ever seen a geological structure like that? If not than probably a rethinking about the parameters might be neeeded.

7) "Is there a standard or beautiful shape of sample variogram against lag intervals?  What is your expects when you do the sample variogram? Does the shape need to be smooth? Is there a standard or theory value for this?"

There are 3 models types that are very used in geostatistics (spherical, gaussian and exponential). They're very used because it's very common to appear variogram structures that resemble them but there's no true rule for this. Some phenomena (like precipitation over time) is cyclical and you get a wave like behaviour in the variogram. Most authors transform the mentioned models with a sin wave before proceeding to kriging or simulation. I theory every model can be used, even those that do not follow the rule "the more the distance, the less correlated".

Final comment: Remember geostatistics is mostly a field of statistics. Things like variograms are tools not the purpose itself. Also I find it useful to think of directional variograms for what they are: subsets of the full set. If the subset is to small than probably the variogram points will have a different variance (let's make the sill) of the global one and suddenly you start to get points that never go above sill (because the variance of your sub-set is less than the one from the global set). If your subset is to big than your start to loose specific information (because you're no longer being selective). For this reason I usually start with small subsets (more lags, less tolerance and so on...) and once I've mapped correctly the variogram in all 180 degress I start to make them bigger (more representative) and closer to the global set, because I still need to find my global anisotropy model.

Hope it helps.

Dear Pedro Correia

Thanks very much for the methodical details and piratical experience about the discussion. I learned a lot from your experience and I believe that your answer and summary will  also conducive to  many other readers. 

Hi Bingwei,

First of all, your questions are very important for optimization geostatistical analysis results.

For the number Lag bins, there is no specific rule to defined the lag size and number of lags. I think the good way is to fit the semivariogram at small lags (distances), then use different bin sizes and refit semivariogram model by changing the lag size and number of lags.

For the search radius of search neighborhood, It is a common practice to limit the data used by defining a circle (or ellipse) to enclose the points that are used to predict values at unmeasured locations. Additionally, to avoid bias in a particular direction, the circle (or ellipse) can be divided into sectors from which an equal number of points are selected. you can specify the number of points, the radius (major/minor axis), and the number of sectors of the circle (or ellipse) to be used for prediction.

I think the next reference will be helpful for the rest questions:

A Practical Guide to Geostatistical Mapping

http://spatial-analyst.net/book/download?size=thumbnail

In addition, with respect to anisotropic variograms, I recommend considering the physical setting when thinking about fitting directional variograms.  For example regardless of the number of lag pairs available, in a riverine setting the only variograms of interest are the long-flow and cross-flow directions.  And based on my experience with many sediment data sets from large navigable rivers the range in the cross flow direction is usually about half the channel width and the long flow variogram cannot be generalized and must be estimated from data. 

The point is that the physical/geological  setting may provide strong clues as to the nature of the variogram that should be most probably.  

Thank you very much for all your answers, here I summarized a workflow for Kriging (attached file). Any comments and suggestions are welcome. 

Thank you in advance!

Hi,

Nice flowchart Bingwei. Local anisotropies procedures may be of great help when doing complex case studies. John Kerns mentioned considering the physical set as a guide for choosing the best models but it is actually possible to use the physical set itself if you have a map of it (and after some processing).

Typically in kriging you use only one variogram ellipsoid for all the mesh (so just one point to point, and/or point to node covariance/correlogram matrix). But with local anisotropies you give a specific ellipsoid for every node (I would advise using stochastic simulation with this procedure to better understand its mathematical uncertainty).

This article: http://link.springer.com/chapter/10.1007/978-90-481-2322-3_15

Gives you such an experiment. If memory does not fail me I believe they've used a map from a fluid simulation to build their local anisotropy maps. So for every node there's the three angles (Azimuth, Dip and Rake) and ratios (main/minor1,main/minor2). Using this they were able to estimate the occurrence of a contaminant by using real samples and the phenomena they knew that caused it's dispersion.

I've also seen somewhere these procedures being used for seismic inversion but I can't really find the article name right now. Nevertheless there's a reason why this algorithms aren't used more often. They are very computer heavy. Most geostatistics algorithms increase their performance by building all covariance values at priori. With this it's no longer possible since the ellipsoid changes every node.

Hi,Pedro, Happy new year! At bellowing question link, professor Donald Myers had given a answer about local anisotropy or non-geometric anisotropy. He mentioned that " Essentially any of the commonly used geostatistics packages incorporate geometric anisotropy but none include non-geometric anisotropy. " For the local anisotropy, you aslo mentioned that, It is heavy computing due to changing correlogram. So we have to try stochastic simulation for improve accuracy. Have I got that right? Or are there any software really can do this?

 https://www.researchgate.net/post/How_can_I_calculate_and_show_the_3D_variogram_surface_to_sub-land-surface_temperature

Happy new year to you too!!!

To my knowledge what professor Donald told you is correct. I also don't know of any geostatistics package that has non-geometric anisotropy procedures. The one I've used was developed within my research center. Let me check (in a couple of days) if I can find you a workable version of that algorithm (the one used in the article I've mentioned).

I've mentioned stochastic simulation not because it's more accurate (as a result), but because it characterizes uncertainty in a far more complete fashion. For instance if you calculate a kriging variance map for, let's say ordinary kriging, you'll see that the lower variances are typically closer to the hard-data. But once you start to insert other criteria in the kriging procedure (such as local anisotropy, or collocated information, or local correlation coefficients, etc.) the variance map will be far more complex. With stochastic simulation instead of having just a mean (kriging mean) and variance (kriging variance) for each node you'll have a distribution (also for each node). Evidently from that distribution you'll be able to calculate the mean and variance (which should be closer to kriging).

This may not sound as great advantage when you look at it at first but consider the implications. For each node you have a distribution of values you consider possible (with your own estimation criteria) to be at that location. That means you can use any of those values to optimize for a final result. A typical solution of this kind is stochastic seismic inversion where you simulate possible images and you compare them to the real seismic (there are a few more steps but I'll skip them for the sake of brevity). In some areas the correlation is good, in others bad. So you use the best areas to generate even more possibilities (collocated simulation) until you get something that having all the criteria you've input in simulation is also strongly correlated with the real seismic.

About computer performance the reason why it's heavy is simple. Usually when you implement geostatistics algorithms the first step on the algorithm (after loading all data) is to build a covariance (usually correlogram because it's safer numerically) matrix that gives you information such as: if distance in x is 3 and distance in y is 9 than the covariance value is cov_matrix[3,9] (again there are few more safeguards but this is the gist of it). We can do this in the beginning because our ellipsoid will be the same every node so every time we get distance X =3, and distance y= 9 we can (instead of re-calculating) retrieve the previously calculated number cov_matrix[3,9]. With local anisotropies this is no longer possible because distance X=3 and distance Y=9 may mean any number of possible ellipsoids. The consequence is we have to calculate the covariances for every single node even if the distances are exactly the same as the previous one.

I think, you can get some important information in the following paper.

A tutorial guide to geostatistics: Computing and modelling variograms and kriging, Catena (Impact Factor: 2.48). 02/2014; 113:56–69. DOI: 10.1016/j.catena.2013.09.006 

Link:

https://www.researchgate.net/publication/259126986_A_tutorial_guide_to_geostatistics_Computing_and_modelling_variograms_and_kriging

Article A tutorial guide to geostatistics: Computing and modelling v...

La técnica de muestreo para variogramas uno la puede asumir de dos maneras:

La primera es utilizando el precepto aleatorio puro en donde uno establece una grilla de puntos aleatoria estratificada del sector, con un numero máximo de puntos.

O haciendo un muestreo por sectores de interes, ¿en que consiste? facil uno establece una mascara la cual deje libre el sector de interes y posteriormente estableciendo sobre este un sistema de grilla aleatoria estratificada.

Mira en si la geoestadistica es un arma de dos filos, como depende mucho de la estadística, tu puedes interevenir sobre sus resultados en funcion de la forma como establezcas el muestreo, por lo cual la forma de establecer este, tiene que estar diseñada de tal forma, que muestre lo verídico del fenómeno (por lo cual depende mucho de tu calidad como investigador) saludos....  

Hi Bingwei,

Sorry for taking so long to put here a version of DSS with local anisotropies. I've put it in the datasets on this link:

https://www.researchgate.net/publication/270902858_DSS_with_local_anisotropies_and_test_data

Notice it's an old test version of this procedure and I haven't done but should anyone require more information just ask and I'll try to answer to the best of my ability.

Again sorry for the delay. Kind regards.

Data DSS with local anisotropies and test data

Hi, Pedro, 

Thank you very much for preparing the software package and documenting it.

I am really appreciated about that. I will try your data and software next week.

Thanks again!

Best Regards!

I really appreciate Pedro's effort to answer these questions which are by any standards make up a nice and valuable manuscript. Bingwei, I would encourage you to push forward with this discussion as a manuscript with Pedro as co-author.

Thanks.

Oh! Thanks for the ego boost:) but, please, there's no need. I'm all for open-knowledge, open-source, open-societies, etc. (among other opens). 

It's a pleasure to change info here in research gate with fellow scientists and I end up learning a lot just because some of the questions stimulate me into thinking more about a subject. It's is own reward really. Nevertheless thanks again for the compliment.

Best regards,

Pedro

Thanks for inspiration from Sarah Alam.  Same as Pedro, I am also always facing to open source and reproducible research and always thinking  that how we can contribute more and more to open world and make it more attractive.

Perhaps making a manuscript from a open discussion and publishing it to a open journal is a way we can make open better.

I really would like to see such kind of open cooperation model. This might be a good idea?  Maybe we can do it and create such open draft way.

How do you think， Pedro?

Hi Bingwei,

Sure thing. Seems like a good initiative and you can count me in. What do you suggest?

Hi, Pedro, Thanks. Let's  start this manuscript project.

First of all I will review the background and summarize your discussions, then we could have deeper discussions on the structures of manuscript. I will contact you.

Hi Bingwei,

I have researched some of these questions a bit and found this paper by Webster and Oliver to be helpful. (Webster, R., Oliver, M.A., 1992. Sample adequately to estimate variograms of soil properties.J. Soil Sci. 43 (1), 177–192.) They recommend at least 50-100 points for a reasonable variogram, but point out that it may vary depending on the study area.

I can't really add anything regarding your other questions.

Have you tried the 'automap' package in R? It allows automatic fitting of the variogram. I've only started trying it out...

Cheers,

Matt

Pedro has made a number of good points,I suggest that it is always important to do some exploratory analysis of your data beginning with histograms (maybe different choices on the bin widths), coded plots of the data locations (coded by the data values),  Compute and plot sample variograms with different choices of lag widths and numbers of lags. Then look at cross validation results with different variogram models/variogram parameters, search neighborhood parameters. There are at least six different cross validation statistics that may be useful.

All these ideas and more are discussed in the tutorial

1991, Myers,D.E.Interpolation and Estimation with Spatially Located Data Chemometrics and Intelligent Laboratory Systems 11, 209-228

It uses a free and very easy to use, interactive package called "geo-EAS", it is old (circa 1988) and runs under DOS but has surprisingly good graphics.  Do a search on Google for "geo-EAS", it was released by the US EPA

1987, A. Warrick and D.E. Myers, Optimization of Sampling Locations for Variogram Calculations Water Resources Research 23, 496-500

1991, Myers,D.E., On Variogram Estimation in Proceedings of the First Inter. Conf. Stat. Comp., Cesme, Turkey, 30 Mar.-2 April 1987, Vol II, American Sciences Press, 261-281

The R package "gstat" is very powerful and easy to use (it includes tutorial materials) but not quite so handy for interactive exploration of the data

Addendum

The sample (empirical, experimental) variogram is not an objective, rather it is a tool to determine/decide on a theoretical variogram. You can't use a sample variogram in the kriging equations. A variogram is a conditionally negative definite function, unfortunately it is not easy in general to determine whether a given functional form satisfies this condition hence in practice one uses a model known to satisfy the condition, e.g. spherical, gaussian, exponential. power.

It is mis-leading to think there is an automatic tool to fit a variogram to a (finite) data set, e.g. automap. There are algorithms that will determine the best fit (best fit using some particular criterion) parameters but they will not choose the model type. A variogram is a characteristic of a random function, however the data is at best only a non-random sample from one realization of the random function. Hence the data do not uniquely characterize the random function.

The answer is NO for all of Pedro's questions.

Since you are at Geovariances you must have access to G. Matheron's many publications (both internal reports and those in journals/conference proceedings). It would be good to review those

I didn't make any questions. Did you mean to tell all my answers are wrong?

Neither I suggested the variogram to be objective (much the opposite, I reinforce the idea that the user is making an interpretation), or that any automatic tool could be perfect to perform such an operation (much by the same reason).

As a side note Geovariances  software do provide tools to make automatic fitting but this was made to broaden the user base to include people that may have no knowledge of geostatistics (thus easing the technical difficulty of the softwares). But the manual interpretation tools are still available. Also there are case studies where the variogram behaves well enough to (most of the times) not need human intervention (besides quality control). Downscaling satellite images is one example (as are others that rely on great amounts of data).

Your questions are not "wrong" but rather for most of them the answer is no, e.g. there is no "best" sample variogram. first of all because there is no single criterion for judging "best" and secondly data sets can vary so much in many different ways. It is not absolutely necessary to even use a sample variogram for kriging. If the variogram model is strictly conditionally negative definite then the system of equations will have a unique solution BUT the results may not be useful. The question of whether the results are useful is not a statistical question but rather has to be answered in the context of the particular application.

There is a big difference between saying that the "sample variogram" doesn't have to be objective and and whether it is an objective. In the one case the word "objective" is an adjective and in the other it is a noun.

Your questions or points pertain to the exploratory analysis you should always do for a data set but how the results of that analysis are used is not absolutely determined.

I know of one example (but which ended up being essentially non conclusive) pertained to studying a potential site for nuclear waste disposal. The site was in an underground old lava flow. A critical question was the vertical thickness of the flow as well as the horizontal thickness. the site was 3000 ft down. Hence the only way to collect data was to drill holes but each such hole was a potential leakage path and drilling the holes was expensive. They also had some outcrop information. the net result was eleven data points which they wanted to use to interpolate the thicknesses, obviously there is some correlation between vertical and horizontal thickness (a further complication).

They used kriging but ultimately decided the results were too unreliable.

Your questions are useful as a guide in the exploratory analysis but the answers are not absolute.

The most suitable Kriging variogram models are Spherical, Circular, and Pentaspherical models.

Hossein

You have not given any crriteria for judging "most suitable" so your answer is meaningless. "most suitable" has to be judged against the usefulness of the results obtained by kriging (and there is not single criterion for that) for a particular data set and the intended use of the results.The variogram model is not the only choice that affects the kriging results.

Donald Myers you are well said.

Thanks a lot

An important step in geostatistical prediction is the selection of suitable correlation structure to model the variogram and its parameters estimation. The predictive performance of each fitted model (spherical, exponential, linear, and Gaussian) can be checked on the basis of cross-validation tests. The values of mean error (ME), mean square error (MSE), reduced mean error (RME), and reduced variance (RV) are estimate to ascertain the performance of the best-fit theoretical models.

Sorry F. Bouteldjaoui but you have not added anything new. Also you have overlooked the possibility that the different statistics you mention may conflict in terms of a conclusion. The mean square error is not a useful statistic.

When you use cross validation you have to make certain choices and the various possible statistics are not equally sensitive to changes in those choices, e.g. parameters for the search neighborhood.

see 1991, Myers,D.E., On Variogram Estimation in Proceedings of the First Inter. Conf. Stat. Comp., Cesme, Turkey, 30 Mar.-2 April 1987, Vol II, American Sciences Press, 261-281