Multidimensional kriging estimation with noisy data observations can be performed using “krig”
[d_est,d_var,lambda_sk,K_dd,k_du,inhood]=krig(pos_known,val_known,pos_est,V,options);
ndata : number of data observations
ndims : dimensions of data location (>=1)
nest : number of data locations to be estimated
pos_known [ndata,ndims] : Locations of data observations
val_known [ndata,1 or 2] : col1 : Data value as measured at 'pos_known'
col2 : Data uncertainty as measured at
'pos_known' (optional)
pos_est [1 ,ndims] : Location of data to be estimated
V : Variogram model, e.g. '1 Sph(100)'
options : kriging options.
Please note data uncertainty is listed as the second column of val_known. Data uncertainty is defined as the variance of Gaussian noise distribution associated to each data measurement. Specified in this manner, the noise on data observations is considered uncorrelated. If only one column is specified, data observations are treated as noise free.
options, the last argument, is a Matlab structure that controls most aspect of the kriging, such as choosing the neighborhood, and kriging method. The names for most of these options are same as the names used to control GSTAT, see “GSTAT from Matlab”
Table 2.1. kriging options
| option.[?] | GSTAT equiv | description | link |
|---|---|---|---|
| sk_mean,mean | sk_mean | Specify the global mean, as used by simple kriging | |
| trend | trend | Krig only the trend. | |
| max | max | The maximum number of data to use in the search neighborhood | |
| polytrend | d | The polynomial order of the trend for each dimension | |
| xvalid | xvalid | If specified as xvalid=1 cross validation on known data locations is performed | |
| isorange | isorange=1 assumes no rotation in anisotropy. If chosen, each entry in the range section of a variogram corresponds to the range for a certain dimension. For example 'Sph(1,10,1,100)' corresponds to a range of 1 in the first dimension, 10 in the second dimension, 1 in the third dimension and 100 in the fourth dimension. If GSTAT is used for kriging in 1, 2 and 3 dimensions, the range selections are translates properly to GSTAT format using the isorange function. |
the Kriging methods: simple Kriging, ordinary Kriging and Kriging with a trend are all available using krig. The difference between these three methods is the way the trend is modelled. Simple Kriging assumed a constant and known mean. Ordinary kriging makes use of an varying unknown mean (that is estimated within the neighborhood). Kriging with a trend model the trend as a smoothly varying polynomial trend. Thus ordinary kriging is but a form of kriging with a trend,; a 0th order polynomial trend.
By default ordinary Kriging in an exhaustive neighborhood (i.e. all data measurements are used all the time) is performed. This is identical to simple Kriging when options.mean=mean(pos_known(:,1))
For each dimension (direction) the order of the polynomial fit of the trend can be specified by options.polyfit
if the value of options.polyfit is an integer scalar, then the same polynomial fit is used for all dimensions. Thus the default linear trend in all directions is identical to specifying options.polyfit=1.
if the value of options.polyfit is an array of integers, each integer in the array must specify the order of the trend for each dimension. The length of the array must be equal to the size of the dimension of the observation. Thus, for 2D data observations one can specify a 2nd order polynomial trend in the first direction and a 0th order fit in the second direction as options.polyfit=[2 0].
By default all data observations are always used. The kriging neighborhood denote the data that are actually used by the kriging systems.
A maximum number of data used by the kriging system is selected as options.max=10. The data locations closes to the point being estimated is retained in the data neighborhood.
See also “nhood” that controls the data neighborhood behavior.
Based on the data observations below, a number of example are shown illustrating the use of kriging in mGstat
In the following example 3 data measurements has been made : x(1)=0; x(5)=3; x(2)=2. Using a Spherical semivariogram model with a range of '0.2' and a sill of '1', the mean and variance of the distribution of the local probability density function at x(2) is found:
pos_known=[1;5;10];
val_known=[0;3;2]; % adding some uncertainty
V='1 Sph(.2)'; % Select variogram model
pos_est=[2]';
[d_est,d_var]=krig(pos_known,val_known,pos_est,V)
d_est =
1.6667
d_var =
1.3333
from: mGstat/examples/mgstat_examples/krig_ex1.m
pos_known=[1;5;10]; val_known=[0;3;2]; % V='1 Sph(.2)'; % Select variogram model pos_est=[0:.1:10]'; [d_est,d_var]=krig(pos_known,val_known,pos_est,V); plot(pos_est,d_est,'k',pos_known,val_known,'ro') print -dpng krigex2
pos_known=[1;5;10]; %
val_known=[0;3;2]; %
V='1 Sph(.2)'; % Select variogram model
pos_est=[0:.1:10]';
[d_est_ok,d_var_ok]=krig(pos_known,val_known,pos_est,V);
options.mean=2;
[d_est_sk,d_var_sk]=krig(pos_known,val_known,pos_est,V,options);
options.trend=1;
[d_est_kt,d_var_kt]=krig(pos_known,val_known,pos_est,V,options);
plot(pos_est,[d_est_sk,d_est_ok,d_est_kt],'-',pos_known,val_known,'ro')
legend('SK','OK','KT','Data')
print -dpng krigex3
rand('seed',1)
ndata=10;
pos_known=rand(ndata,1)*10;
val_known=randn(ndata,1); %
V='1 Sph(.2)'; % Select variogram model
pos_est=[0:.1:10]';
clear options;
[d_est_ok,d_var_ok]=krig(pos_known,val_known,pos_est,V,options);
options.mean=2;
[d_est_sk,d_var_sk]=krig(pos_known,val_known,pos_est,V,options);
options.trend=1;
[d_est_kt,d_var_kt]=krig(pos_known,val_known,pos_est,V,options);
plot(pos_est,[d_est_sk,d_est_ok,d_est_kt],'-',pos_known,val_known,'k*')
legend('SK','OK','KT','Data')
print -dpng krigex4
rand('seed',1)
ndata=10;
pos_known=rand(ndata,1)*10;
val_known=randn(ndata,1); %
V='1 Sph(.2)'; % Select variogram model
pos_est=[0:.1:10]';
clear options;options.max=4;
[d_est_ok,d_var_ok]=krig(pos_known,val_known,pos_est,V,options);
options.mean=2;
[d_est_sk,d_var_sk]=krig(pos_known,val_known,pos_est,V,options);
options.trend=1;
[d_est_kt,d_var_kt]=krig(pos_known,val_known,pos_est,V,options);
plot(pos_est,[d_est_sk,d_est_ok,d_est_kt],'-',pos_known,val_known,'k*')
legend('SK','OK','KT','Data')
print -dpng krigex5
rand('seed',1)
ndata=30;
pos_known=rand(ndata,1)*10;
val_known=randn(ndata,1); %
val_var=zeros(ndata,1)+.1; %
V='1 Sph(.1)'; % Select variogram model
pos_est=[-2:.1:12]';
clear options;options.max=4;
[d_est_ok,d_var_ok]=krig(pos_known,[val_known val_var],pos_est,V,options);
options.mean=mean(val_known);
[d_est_sk,d_var_sk]=krig(pos_known,[val_known val_var],pos_est,V,options);
options.trend=1;
[d_est_kt,d_var_kt]=krig(pos_known,[val_known val_var],pos_est,V,options);
plot(pos_est,[d_est_sk,d_est_ok,d_est_kt],'-',pos_known,val_known,'k*')
legend('SK','OK','KT','Data')
print -dpng krigex6
| © 2004-2010 Thomas Mejer Hansen. | This site is hosted by
|