Standard 2-D case study¶

Import packages¶

In [1]:
import numpy as np
import gstlearn as gl
import gstlearn.plot as gp
import gstlearn.document as gdoc

Global variables

In [2]:
verbose = True
graphics = True
gl.OptCst.define(gl.ECst.NTCOL, 6)
gdoc.setNoScroll()

Reading data¶

The data are stored in a CSV format in the file called Pollution.dat

In [3]:
filepath = gdoc.loadData("Pollution", "Pollution.dat")
mydb = gl.Db.createFromCSV(filepath, gl.CSVformat())
mydb.setLocators(["X", "Y"], gl.ELoc.X)
mydb.setLocator("Zn", gl.ELoc.Z)
if verbose:
    dbfmt = gl.DbStringFormat()
    dbfmt.setFlags(flag_resume=True, flag_vars=True, flag_extend=True)
    mydb.display(dbfmt)
Data Base Characteristics
=========================

Data Base Summary
-----------------
File is organized as a set of isolated points
Space dimension              = 2
Number of Columns            = 5
Total number of samples      = 102

Data Base Extension
-------------------
Coor #1 - Min =     109.850 - Max =     143.010 - Ext = 33.16
Coor #2 - Min =     483.660 - Max =     513.040 - Ext = 29.38

Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = X - Locator = x1
Column = 2 - Name = Y - Locator = x2
Column = 3 - Name = Zn - Locator = z1
Column = 4 - Name = Pb - Locator = p1

Accessing to the variable names

In [4]:
print("List of all variable names =", mydb.getAllNames())
List of all variable names = ('rank', 'X', 'Y', 'Zn', 'Pb')

Extracting the vector containing the Zn variable in order to perform a selection

In [5]:
tabZn = mydb.getColumn("Zn")  # equivalent to mydb["Zn"] or mydb[3]
selZn = tabZn < 20
mydb.addSelection(selZn, "sel")
mydb.setLocator("Pb", gl.ELoc.Z)
dbfmt = gl.DbStringFormat.createFromFlags(flag_stats=True)
if verbose:
    mydb.display(dbfmt)
Data Base Characteristics
=========================

Data Base Summary
-----------------
File is organized as a set of isolated points
Space dimension              = 2
Number of Columns            = 6
Total number of samples      = 102
Number of active samples     = 99

Data Base Statistics
--------------------
1 - Name rank - Locator NA
 Nb of data          =         102
 Nb of active values =          99
 Minimum value       =       1.000
 Maximum value       =     102.000
 Mean value          =      51.808
 Standard Deviation  =      29.114
 Variance            =     847.650
2 - Name X - Locator x1
 Nb of data          =         102
 Nb of active values =          99
 Minimum value       =     109.850
 Maximum value       =     143.010
 Mean value          =     119.924
 Standard Deviation  =       6.549
 Variance            =      42.890
3 - Name Y - Locator x2
 Nb of data          =         102
 Nb of active values =          99
 Minimum value       =     483.660
 Maximum value       =     513.040
 Mean value          =     498.622
 Standard Deviation  =       8.000
 Variance            =      64.001
4 - Name Zn - Locator NA
 Nb of data          =         102
 Nb of active values =          99
 Minimum value       =       1.090
 Maximum value       =      12.100
 Mean value          =       2.881
 Standard Deviation  =       1.667
 Variance            =       2.779
5 - Name Pb - Locator z1
 Nb of data          =         102
 Nb of active values =          99
 Minimum value       =       3.000
 Maximum value       =      12.700
 Mean value          =       5.658
 Standard Deviation  =       1.697
 Variance            =       2.881
6 - Name sel - Locator sel
 Nb of data          =         102
 Nb of active values =          99
 Minimum value       =       1.000
 Maximum value       =       1.000
 Mean value          =       1.000
 Standard Deviation  =       0.000
 Variance            =       0.000

Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = X - Locator = x1
Column = 2 - Name = Y - Locator = x2
Column = 3 - Name = Zn - Locator = NA
Column = 4 - Name = Pb - Locator = z1
Column = 5 - Name = sel - Locator = sel

Display my Data (with samples represented by color and size)

In [6]:
if graphics:
    gp.plot(mydb, nameColor="Pb")
    gp.decoration(title="Data Set")
No description has been provided for this image

Variograms¶

We first define the geometry of the variogram calculations

In [7]:
myVarioParamOmni = gl.VarioParam()
mydir = gl.DirParam.create(nlag=10, dlag=1.0)
myVarioParamOmni.addDir(mydir)

We use the variogram definition in order to calculate the variogram cloud.

In [8]:
dbcloud = gl.db_vcloud(mydb, myVarioParamOmni)

We recall that the Variogram cloud is calculated by filling an underlying grid where each cell is painted according to the number of pairs at the given distance and given variability. Representing the variogram cloud

In [9]:
if graphics:
    gp.plot(dbcloud, "Cloud*")
    gp.decoration(title="Variogram Cloud")
No description has been provided for this image

Calculating the experimental omni-directional variogram

In [10]:
myVarioOmni = gl.Vario(myVarioParamOmni)
err = myVarioOmni.compute(mydb, gl.ECalcVario.VARIOGRAM)
if verbose:
    myVarioOmni.display()
Variogram characteristics
=========================
Number of variable(s)       = 1
Number of direction(s)      = 1
Space dimension             = 2
Variable(s)                 = [Pb]

Variance-Covariance Matrix      2.881

Direction #1
------------
Number of lags              = 10
Direction coefficients      =       1.000      0.000
Direction angles (degrees)  =       0.000
Tolerance on direction      =      90.000 (degrees)
Calculation lag             =       1.000
Tolerance on distance       =      50.000 (Percent of the lag value)

For variable 1
       Rank     Npairs   Distance      Value
          0      3.000      0.389      0.462
          1    123.000      1.081      1.495
          2    183.000      2.038      1.620
          3    205.000      3.006      2.526
          4    231.000      4.013      2.240
          5    229.000      5.036      2.524
          6    198.000      5.962      2.396
          7    187.000      7.000      2.708
          8    204.000      7.996      2.772
          9    184.000      8.990      2.868

The variogram is represented graphically for a quick check

In [11]:
if graphics:
    gp.plot(myVarioOmni)
    gp.decoration(title="Omni-directional Variogram for Pb")
No description has been provided for this image

Calculate a variogram in several directions

In [12]:
myvarioParam = gl.VarioParam()
mydirs = gl.DirParam.createMultiple(4, 10, 1.0)
myvarioParam.addMultiDirs(mydirs)
myvario = gl.Vario(myvarioParam)
myvario.compute(mydb, gl.ECalcVario.VARIOGRAM)
if verbose:
    myvario.display()
Variogram characteristics
=========================
Number of variable(s)       = 1
Number of direction(s)      = 4
Space dimension             = 2
Variable(s)                 = [Pb]

Variance-Covariance Matrix      2.881

Direction #1
------------
Number of lags              = 10
Direction coefficients      =       1.000      0.000
Direction angles (degrees)  =       0.000
Tolerance on direction      =      22.500 (degrees)
Calculation lag             =       1.000
Tolerance on distance       =      50.000 (Percent of the lag value)

For variable 1
       Rank     Npairs   Distance      Value
          0      1.000      0.410      0.180
          1     29.000      1.094      1.634
          2     47.000      2.079      1.415
          3     53.000      3.003      2.824
          4     63.000      3.999      2.348
          5     66.000      5.035      2.319
          6     60.000      5.978      3.115
          7     52.000      7.045      2.746
          8     52.000      8.020      3.927
          9     37.000      8.980      2.554

Direction #2
------------
Number of lags              = 10
Direction coefficients      =       0.707      0.707
Direction angles (degrees)  =      45.000
Tolerance on direction      =      22.500 (degrees)
Calculation lag             =       1.000
Tolerance on distance       =      50.000 (Percent of the lag value)

For variable 1
       Rank     Npairs   Distance      Value
          0      1.000      0.344      0.080
          1     31.000      1.051      1.113
          2     50.000      1.960      1.890
          3     62.000      2.999      2.443
          4     58.000      4.014      2.701
          5     51.000      5.016      2.702
          6     36.000      5.999      1.833
          7     37.000      7.015      2.130
          8     50.000      7.997      2.060
          9     53.000      8.995      2.381

Direction #3
------------
Number of lags              = 10
Direction coefficients      =       0.000      1.000
Direction angles (degrees)  =      90.000
Tolerance on direction      =      22.500 (degrees)
Calculation lag             =       1.000
Tolerance on distance       =      50.000 (Percent of the lag value)

For variable 1
       Rank     Npairs   Distance      Value
          1     32.000      1.149      1.631
          2     39.000      2.080      1.670
          3     39.000      2.979      2.511
          4     48.000      4.012      2.120
          5     51.000      5.029      3.055
          6     47.000      5.939      2.856
          7     49.000      6.965      2.386
          8     42.000      7.952      2.708
          9     41.000      9.018      2.320

Direction #4
------------
Number of lags              = 10
Direction coefficients      =      -0.707      0.707
Direction angles (degrees)  =     135.000
Tolerance on direction      =      22.500 (degrees)
Calculation lag             =       1.000
Tolerance on distance       =      50.000 (Percent of the lag value)

For variable 1
       Rank     Npairs   Distance      Value
          0      1.000      0.411      1.125
          1     31.000      1.028      1.606
          2     47.000      2.044      1.496
          3     51.000      3.040      2.330
          4     62.000      4.028      1.791
          5     61.000      5.058      2.155
          6     55.000      5.939      1.587
          7     49.000      6.975      3.425
          8     60.000      8.004      2.408
          9     53.000      8.972      3.996
In [13]:
if graphics:
    gp.plot(myvario, idir=-1)
    gp.decoration(title="Multi-Directional Variogram of Pb")
No description has been provided for this image

Calculating the Variogram Map

In [14]:
myvmap = gl.db_vmap(mydb, calculType=gl.ECalcVario.VARIOGRAM, nxx=[20, 20])
if verbose:
    myvmap.display()
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 5
Total number of samples      = 1681

Grid characteristics:
---------------------
Origin :     -33.160    -29.380
Mesh   :       1.658      1.469
Number :          41         41

Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = x1 - Locator = x1
Column = 2 - Name = x2 - Locator = x2
Column = 3 - Name = VMAP.Pb.Var - Locator = z1
Column = 4 - Name = VMAP.Pb.Nb - Locator = NA
In [15]:
if graphics:
    gp.plot(myvmap, "*Var")
    gp.decoration(title="Variogram Map")
No description has been provided for this image

Model¶

Fitting a Model. We call the Automatic Fitting procedure providing the list of covariance functions to be tested.

In [16]:
mymodel = gl.Model.createFromDb(mydb)
err = mymodel.fit(myvario, [gl.ECov.EXPONENTIAL, gl.ECov.SPHERICAL])

Visualizing the resulting model, overlaid on the experimental variogram

In [17]:
if graphics:
    gp.varmod(myvario, mymodel)
    gp.decoration(title="Model for Pb")
No description has been provided for this image

Model with equality constraints¶

We can impose some constraints on the parameters during the fit. For instance here, we impose an equality constraint on the range (range = 1).

In [18]:
myModelConstrained = gl.Model.createFromDb(mydb)
constr = gl.Constraints()
paramid = gl.CovParamId(0, 0, gl.EConsElem.RANGE, 0, 0)
constr.addItem(gl.ConsItem(paramid, gl.EConsType.EQUAL, 1.0))
err = myModelConstrained.fit(myVarioOmni, [gl.ECov.SPHERICAL], constr)
if err > 0:
    print("Error while fitting the model")
myModelConstrained
Out[18]:
Model characteristics
=====================
Space dimension              = 2
Number of variable(s)        = 1
Number of basic structure(s) = 1
Number of drift function(s)  = 0
Number of drift equation(s)  = 0

Covariance Part
---------------
Spherical
- Sill         =       2.101
- Range        =       1.000
Total Sill     =       2.101
Known Mean(s)      0.000

We can impose inequality constraints by using EConsType.LOWER or EConsType.UPPER.

Adding a drift¶

In [19]:
mymodel.setDriftIRF()
if verbose:
    mymodel.display()
Model characteristics
=====================
Space dimension              = 2
Number of variable(s)        = 1
Number of basic structure(s) = 2
Number of drift function(s)  = 1
Number of drift equation(s)  = 1

Covariance Part
---------------
Exponential
- Sill         =       1.039
- Ranges       =       2.103      0.386
- Theo. Ranges =       0.702      0.129
- Angles       =      44.853      0.000
- Rotation Matrix
                 [,  0]     [,  1]
      [  0,]     -0.709      0.705
      [  1,]     -0.705     -0.709
Spherical
- Sill         =       1.606
- Ranges       =       6.909      5.009
- Angles       =     136.493      0.000
- Rotation Matrix
                 [,  0]     [,  1]
      [  0,]     -0.725     -0.688
      [  1,]      0.688     -0.725
Total Sill     =       2.645

Drift Part
----------
Universality_Condition

Defining the Neighborhood¶

We initiate a Neigborhood (Moving with a small number of samples for Demonstration)

In [20]:
myneigh = gl.NeighMoving.create(flag_xvalid=False, nmaxi=6, radius=10)
if verbose:
    myneigh.display()
Moving Neighborhood
===================
Minimum number of samples           = 1
Maximum number of samples           = 6
Maximum horizontal distance         = 10

Checking the Moving Neighborhood¶

We must first create a Grid which covers the area of interest

In [21]:
mygrid = gl.DbGrid.createCoveringDb(mydb, [], [0.5, 0.5], [], [2, 2])
if verbose:
    mygrid.display()
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 2
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2

We can now test the neighborhood characteristics for each node of the previously defined grid.

In [22]:
err = gl.test_neigh(mydb, mygrid, mymodel, myneigh)
if err > 0:
    print("Error while running test_neigh")
if verbose:
    mygrid.display()
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 7
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = z1

We can visualize some of the newly created variables, such as:

  • the number of points per neighborhood
In [23]:
if graphics:
    gp.plot(mygrid, "Neigh*Number")
    gp.decoration(title="Number of Samples per Neighborhood")
No description has been provided for this image
  • the one giving the maximum distance per neighborhood
In [24]:
if graphics:
    gp.plot(mygrid, "Neigh*MaxDist")
    gp.decoration(title="Maximum Distance per Neighborhood")
No description has been provided for this image

Cross-validation¶

We can now process the cross-validation step

In [25]:
err = gl.xvalid(mydb, mymodel, myneigh)
if err > 0:
    print("Error while running xvalid")
if verbose:
    mydb.display()
Data Base Characteristics
=========================

Data Base Summary
-----------------
File is organized as a set of isolated points
Space dimension              = 2
Number of Columns            = 8
Total number of samples      = 102
Number of active samples     = 99

Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = X - Locator = x1
Column = 2 - Name = Y - Locator = x2
Column = 3 - Name = Zn - Locator = NA
Column = 4 - Name = Pb - Locator = NA
Column = 5 - Name = sel - Locator = sel
Column = 6 - Name = Xvalid.Pb.esterr - Locator = z1
Column = 7 - Name = Xvalid.Pb.stderr - Locator = NA
In [26]:
if graphics:
    gp.histogram(mydb, name="Xvalid.Pb.stderr", bins=50, range=(-4, 3))
No description has been provided for this image

Estimation by Kriging¶

We now perform the Estimation by Ordinary Kriging. The Neighborhood is changed into a Unique Neighborhood.

In [27]:
mydb.setLocator("Pb", gl.ELoc.Z)
myneigh = gl.NeighUnique.create()
err = gl.kriging(mydb, mygrid, mymodel, myneigh)
if err > 0:
    print("Error while running kriging")
if verbose:
    mygrid.display(dbfmt)
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 9
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Data Base Statistics
--------------------
1 - Name x1 - Locator x1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     107.850
 Maximum value       =     144.850
 Mean value          =     126.350
 Standard Deviation  =      10.824
 Variance            =     117.167
2 - Name x2 - Locator x2
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     481.660
 Maximum value       =     515.160
 Mean value          =     498.410
 Standard Deviation  =       9.814
 Variance            =      96.313
3 - Name Neigh.Pb.Number - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       6.000
 Mean value          =       5.381
 Standard Deviation  =       1.538
 Variance            =       2.366
4 - Name Neigh.Pb.MaxDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.935
 Maximum value       =       9.999
 Mean value          =       5.598
 Standard Deviation  =       2.497
 Variance            =       6.233
5 - Name Neigh.Pb.MinDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.014
 Maximum value       =       9.978
 Mean value          =       3.484
 Standard Deviation  =       2.514
 Variance            =       6.320
6 - Name Neigh.Pb.NbNESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       1.000
 Mean value          =       1.000
 Standard Deviation  =       0.000
 Variance            =       0.000
7 - Name Neigh.Pb.NbCESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.000
 Maximum value       =       0.000
 Mean value          =       0.000
 Standard Deviation  =       0.000
 Variance            =       0.000
8 - Name Kriging.Pb.estim - Locator z1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.408
 Maximum value       =      11.493
 Mean value          =       6.114
 Standard Deviation  =       0.641
 Variance            =       0.410
9 - Name Kriging.Pb.stdev - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       0.234
 Maximum value       =       1.655
 Mean value          =       1.540
 Standard Deviation  =       0.160
 Variance            =       0.026

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = NA
Column = 7 - Name = Kriging.Pb.estim - Locator = z1
Column = 8 - Name = Kriging.Pb.stdev - Locator = NA

Visualizing the results

In [28]:
if graphics:
    gp.plot(mygrid, "Kriging.Pb.estim")
    gp.plot(mydb, nameColor="Pb")
    gp.decoration(title="Estimate of Pb")
No description has been provided for this image
In [29]:
if graphics:
    gp.plot(mygrid, "Kriging.Pb.stdev")
    gp.plot(mydb, nameColor="Pb")
    gp.decoration(title="St. Deviation of Pb")
No description has been provided for this image

Simulations¶

We must first transform the Data into Gaussian

In [30]:
myanamPb = gl.AnamHermite(30)
myanamPb.fitFromLocator(mydb)
if verbose:
    myanamPb.display()
Hermitian Anamorphosis
----------------------
Minimum absolute value for Y  = -2.7
Maximum absolute value for Y  = 2.6
Minimum absolute value for Z  = 3.0029
Maximum absolute value for Z  = 12.9777
Minimum practical value for Y = -2.7
Maximum practical value for Y = 2.6
Minimum practical value for Z = 3.0029
Maximum practical value for Z = 12.9777
Mean                          = 5.65758
Variance                      = 2.86296
Number of Hermite polynomials = 30
Normalized coefficients for Hermite polynomials (punctual variable)
                 [,  0]     [,  1]     [,  2]     [,  3]     [,  4]     [,  5]     [,  6]
      [  0,]      5.658     -1.625      0.440     -0.069     -0.017      0.082     -0.061
      [  7,]      0.001      0.036     -0.044      0.004      0.047     -0.030     -0.029
      [ 14,]      0.037      0.007     -0.031      0.010      0.018     -0.019     -0.003
      [ 21,]      0.019     -0.010     -0.014      0.019      0.006     -0.023      0.004
      [ 28,]      0.022     -0.013

We can produce the Gaussian Anamorphosis graphically within its definition domain.

In [31]:
if graphics:
    gp.plot(myanamPb)
No description has been provided for this image

The next step consists in translating the target variable ('Pb') into its Gaussian transform

In [32]:
mydb.setLocator("Pb", gl.ELoc.Z)
err = myanamPb.rawToGaussianByLocator(mydb)
if verbose:
    mydb.display()
Data Base Characteristics
=========================

Data Base Summary
-----------------
File is organized as a set of isolated points
Space dimension              = 2
Number of Columns            = 9
Total number of samples      = 102
Number of active samples     = 99

Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = X - Locator = x1
Column = 2 - Name = Y - Locator = x2
Column = 3 - Name = Zn - Locator = NA
Column = 4 - Name = Pb - Locator = NA
Column = 5 - Name = sel - Locator = sel
Column = 6 - Name = Xvalid.Pb.esterr - Locator = NA
Column = 7 - Name = Xvalid.Pb.stderr - Locator = NA
Column = 8 - Name = Y.Pb - Locator = z1

We quickly calculate experimental (omni-directional) variograms of the gaussian variable (locater Z) using the already defined directions

In [33]:
myvarioParam = gl.VarioParam()
mydir = gl.DirParam(10, 1.0)
myvarioParam.addDir(mydir)
myVario = gl.Vario(myvarioParam)
err = myvario.compute(mydb, gl.ECalcVario.VARIOGRAM)

We fit the model by automatic fit (with the constraints that the total sill be equal to 1).

In [34]:
mymodelG = gl.Model.createFromDb(mydb)
err = mymodelG.fit(myvario, [gl.ECov.EXPONENTIAL])
if graphics:
    ax = gp.varmod(myvario, mymodelG)
    ax.decoration(title="Model for Gaussian Pb")
No description has been provided for this image

We perform a set of 10 conditional simulations using the Turning Bands Method.

In [35]:
err = gl.simtub(mydb, mygrid, mymodel, myneigh, nbsimu=10)
if verbose:
    mygrid.display()
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 19
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = NA
Column = 7 - Name = Kriging.Pb.estim - Locator = NA
Column = 8 - Name = Kriging.Pb.stdev - Locator = NA
Column = 9 - Name = Simu.Y.Pb.S1 - Locator = z1
Column = 10 - Name = Simu.Y.Pb.S2 - Locator = z2
Column = 11 - Name = Simu.Y.Pb.S3 - Locator = z3
Column = 12 - Name = Simu.Y.Pb.S4 - Locator = z4
Column = 13 - Name = Simu.Y.Pb.S5 - Locator = z5
Column = 14 - Name = Simu.Y.Pb.S6 - Locator = z6
Column = 15 - Name = Simu.Y.Pb.S7 - Locator = z7
Column = 16 - Name = Simu.Y.Pb.S8 - Locator = z8
Column = 17 - Name = Simu.Y.Pb.S9 - Locator = z9
Column = 18 - Name = Simu.Y.Pb.S10 - Locator = z10

Some statistics on the Conditional simulations in Gaussian scale

In [36]:
mygrid.deleteColumn("Stats.*")
mygrid.statisticsBySample(
    ["Simu.Y.*"],
    [
        gl.EStatOption.MINI,
        gl.EStatOption.MAXI,
        gl.EStatOption.MEAN,
        gl.EStatOption.STDV,
    ],
    True,
)
if verbose:
    dbfmt = gl.DbStringFormat()
    dbfmt.setFlags(flag_resume=True, flag_stats=True)
    mygrid.display(dbfmt)
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 23
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Data Base Statistics
--------------------
1 - Name x1 - Locator x1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     107.850
 Maximum value       =     144.850
 Mean value          =     126.350
 Standard Deviation  =      10.824
 Variance            =     117.167
2 - Name x2 - Locator x2
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     481.660
 Maximum value       =     515.160
 Mean value          =     498.410
 Standard Deviation  =       9.814
 Variance            =      96.313
3 - Name Neigh.Pb.Number - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       6.000
 Mean value          =       5.381
 Standard Deviation  =       1.538
 Variance            =       2.366
4 - Name Neigh.Pb.MaxDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.935
 Maximum value       =       9.999
 Mean value          =       5.598
 Standard Deviation  =       2.497
 Variance            =       6.233
5 - Name Neigh.Pb.MinDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.014
 Maximum value       =       9.978
 Mean value          =       3.484
 Standard Deviation  =       2.514
 Variance            =       6.320
6 - Name Neigh.Pb.NbNESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       1.000
 Mean value          =       1.000
 Standard Deviation  =       0.000
 Variance            =       0.000
7 - Name Neigh.Pb.NbCESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.000
 Maximum value       =       0.000
 Mean value          =       0.000
 Standard Deviation  =       0.000
 Variance            =       0.000
8 - Name Kriging.Pb.estim - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.408
 Maximum value       =      11.493
 Mean value          =       6.114
 Standard Deviation  =       0.641
 Variance            =       0.410
9 - Name Kriging.Pb.stdev - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       0.234
 Maximum value       =       1.655
 Mean value          =       1.540
 Standard Deviation  =       0.160
 Variance            =       0.026
10 - Name Simu.Y.Pb.S1 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -5.146
 Maximum value       =       5.909
 Mean value          =       0.295
 Standard Deviation  =       1.646
 Variance            =       2.708
11 - Name Simu.Y.Pb.S2 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -5.775
 Maximum value       =       5.453
 Mean value          =      -0.094
 Standard Deviation  =       1.546
 Variance            =       2.392
12 - Name Simu.Y.Pb.S3 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -5.496
 Maximum value       =       6.225
 Mean value          =       0.119
 Standard Deviation  =       1.584
 Variance            =       2.509
13 - Name Simu.Y.Pb.S4 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -4.897
 Maximum value       =       6.261
 Mean value          =       0.433
 Standard Deviation  =       1.587
 Variance            =       2.518
14 - Name Simu.Y.Pb.S5 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -4.759
 Maximum value       =       6.695
 Mean value          =       0.420
 Standard Deviation  =       1.588
 Variance            =       2.521
15 - Name Simu.Y.Pb.S6 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -5.215
 Maximum value       =       7.209
 Mean value          =       0.304
 Standard Deviation  =       1.603
 Variance            =       2.569
16 - Name Simu.Y.Pb.S7 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -4.939
 Maximum value       =       6.154
 Mean value          =       0.370
 Standard Deviation  =       1.535
 Variance            =       2.358
17 - Name Simu.Y.Pb.S8 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -5.058
 Maximum value       =       6.995
 Mean value          =       0.323
 Standard Deviation  =       1.568
 Variance            =       2.459
18 - Name Simu.Y.Pb.S9 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -4.184
 Maximum value       =       5.380
 Mean value          =       0.399
 Standard Deviation  =       1.507
 Variance            =       2.271
19 - Name Simu.Y.Pb.S10 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -5.266
 Maximum value       =       6.048
 Mean value          =       0.011
 Standard Deviation  =       1.591
 Variance            =       2.533
20 - Name Stats.MINI - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -5.775
 Maximum value       =       0.847
 Mean value          =      -2.099
 Standard Deviation  =       0.936
 Variance            =       0.876
21 - Name Stats.MAXI - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -1.201
 Maximum value       =       7.209
 Mean value          =       2.681
 Standard Deviation  =       1.024
 Variance            =       1.048
22 - Name Stats.MEAN - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =      -2.125
 Maximum value       =       2.241
 Mean value          =       0.258
 Standard Deviation  =       0.577
 Variance            =       0.333
23 - Name Stats.STDV - Locator z1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       0.267
 Maximum value       =       2.773
 Mean value          =       1.431
 Standard Deviation  =       0.365
 Variance            =       0.134

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = NA
Column = 7 - Name = Kriging.Pb.estim - Locator = NA
Column = 8 - Name = Kriging.Pb.stdev - Locator = NA
Column = 9 - Name = Simu.Y.Pb.S1 - Locator = NA
Column = 10 - Name = Simu.Y.Pb.S2 - Locator = NA
Column = 11 - Name = Simu.Y.Pb.S3 - Locator = NA
Column = 12 - Name = Simu.Y.Pb.S4 - Locator = NA
Column = 13 - Name = Simu.Y.Pb.S5 - Locator = NA
Column = 14 - Name = Simu.Y.Pb.S6 - Locator = NA
Column = 15 - Name = Simu.Y.Pb.S7 - Locator = NA
Column = 16 - Name = Simu.Y.Pb.S8 - Locator = NA
Column = 17 - Name = Simu.Y.Pb.S9 - Locator = NA
Column = 18 - Name = Simu.Y.Pb.S10 - Locator = NA
Column = 19 - Name = Stats.MINI - Locator = NA
Column = 20 - Name = Stats.MAXI - Locator = NA
Column = 21 - Name = Stats.MEAN - Locator = NA
Column = 22 - Name = Stats.STDV - Locator = z1

We visualize a conditional simulation in Gaussian scale

In [37]:
if graphics:
    gp.plot(mygrid, "Simu.Y.Pb.S1")
    gp.plot(mydb, nameColor="Pb")
    gp.decoration(title="One Simulation of Pb in Gaussian Scale")
No description has been provided for this image

We turn the Gaussian conditional simulations into Raw scale (using the Anamorphosis back transform) and get rid of the Gaussian conditional simulations.

In [38]:
myanamPb.gaussianToRaw(mygrid, "Simu.Y.*")
mygrid.deleteColumn("Simu.Y.*")
mygrid.deleteColumn("Stats.*")
if verbose:
    mygrid.display()
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 19
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = NA
Column = 7 - Name = Kriging.Pb.estim - Locator = NA
Column = 8 - Name = Kriging.Pb.stdev - Locator = NA
Column = 9 - Name = Z.Simu.Y.Pb.S1 - Locator = z1
Column = 10 - Name = Z.Simu.Y.Pb.S2 - Locator = z2
Column = 11 - Name = Z.Simu.Y.Pb.S3 - Locator = z3
Column = 12 - Name = Z.Simu.Y.Pb.S4 - Locator = z4
Column = 13 - Name = Z.Simu.Y.Pb.S5 - Locator = z5
Column = 14 - Name = Z.Simu.Y.Pb.S6 - Locator = z6
Column = 15 - Name = Z.Simu.Y.Pb.S7 - Locator = z7
Column = 16 - Name = Z.Simu.Y.Pb.S8 - Locator = z8
Column = 17 - Name = Z.Simu.Y.Pb.S9 - Locator = z9
Column = 18 - Name = Z.Simu.Y.Pb.S10 - Locator = z10

We calculate some statistics on the Conditional Simulations in Raw scale.

In [39]:
mygrid.deleteColumn("Stats.*")
mygrid.statisticsBySample(
    ["Z.Simu.*"],
    [
        gl.EStatOption.MINI,
        gl.EStatOption.MAXI,
        gl.EStatOption.MEAN,
        gl.EStatOption.STDV,
    ],
    True,
)
if verbose:
    dbfmt = gl.DbStringFormat()
    dbfmt.setFlags(flag_resume=True, flag_stats=True)
    mygrid.display(dbfmt)
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 23
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Data Base Statistics
--------------------
1 - Name x1 - Locator x1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     107.850
 Maximum value       =     144.850
 Mean value          =     126.350
 Standard Deviation  =      10.824
 Variance            =     117.167
2 - Name x2 - Locator x2
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     481.660
 Maximum value       =     515.160
 Mean value          =     498.410
 Standard Deviation  =       9.814
 Variance            =      96.313
3 - Name Neigh.Pb.Number - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       6.000
 Mean value          =       5.381
 Standard Deviation  =       1.538
 Variance            =       2.366
4 - Name Neigh.Pb.MaxDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.935
 Maximum value       =       9.999
 Mean value          =       5.598
 Standard Deviation  =       2.497
 Variance            =       6.233
5 - Name Neigh.Pb.MinDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.014
 Maximum value       =       9.978
 Mean value          =       3.484
 Standard Deviation  =       2.514
 Variance            =       6.320
6 - Name Neigh.Pb.NbNESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       1.000
 Mean value          =       1.000
 Standard Deviation  =       0.000
 Variance            =       0.000
7 - Name Neigh.Pb.NbCESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.000
 Maximum value       =       0.000
 Mean value          =       0.000
 Standard Deviation  =       0.000
 Variance            =       0.000
8 - Name Kriging.Pb.estim - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.408
 Maximum value       =      11.493
 Mean value          =       6.114
 Standard Deviation  =       0.641
 Variance            =       0.410
9 - Name Kriging.Pb.stdev - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       0.234
 Maximum value       =       1.655
 Mean value          =       1.540
 Standard Deviation  =       0.160
 Variance            =       0.026
10 - Name Z.Simu.Y.Pb.S1 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.593
 Standard Deviation  =       2.989
 Variance            =       8.931
11 - Name Z.Simu.Y.Pb.S2 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.889
 Standard Deviation  =       2.596
 Variance            =       6.741
12 - Name Z.Simu.Y.Pb.S3 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.252
 Standard Deviation  =       2.771
 Variance            =       7.676
13 - Name Z.Simu.Y.Pb.S4 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.777
 Standard Deviation  =       2.993
 Variance            =       8.957
14 - Name Z.Simu.Y.Pb.S5 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.718
 Standard Deviation  =       2.901
 Variance            =       8.415
15 - Name Z.Simu.Y.Pb.S6 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.600
 Standard Deviation  =       2.911
 Variance            =       8.474
16 - Name Z.Simu.Y.Pb.S7 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.665
 Standard Deviation  =       2.893
 Variance            =       8.371
17 - Name Z.Simu.Y.Pb.S8 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.545
 Standard Deviation  =       2.842
 Variance            =       8.079
18 - Name Z.Simu.Y.Pb.S9 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.700
 Standard Deviation  =       2.870
 Variance            =       8.239
19 - Name Z.Simu.Y.Pb.S10 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.088
 Standard Deviation  =       2.727
 Variance            =       7.439
20 - Name Stats.MINI - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =       6.887
 Mean value          =       3.400
 Standard Deviation  =       0.536
 Variance            =       0.287
21 - Name Stats.MAXI - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.942
 Maximum value       =      12.978
 Mean value          =      11.503
 Standard Deviation  =       2.070
 Variance            =       4.283
22 - Name Stats.MEAN - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.340
 Maximum value       =      10.577
 Mean value          =       6.483
 Standard Deviation  =       1.017
 Variance            =       1.035
23 - Name Stats.STDV - Locator z1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       0.323
 Maximum value       =       4.549
 Mean value          =       2.577
 Standard Deviation  =       0.734
 Variance            =       0.538

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = NA
Column = 7 - Name = Kriging.Pb.estim - Locator = NA
Column = 8 - Name = Kriging.Pb.stdev - Locator = NA
Column = 9 - Name = Z.Simu.Y.Pb.S1 - Locator = NA
Column = 10 - Name = Z.Simu.Y.Pb.S2 - Locator = NA
Column = 11 - Name = Z.Simu.Y.Pb.S3 - Locator = NA
Column = 12 - Name = Z.Simu.Y.Pb.S4 - Locator = NA
Column = 13 - Name = Z.Simu.Y.Pb.S5 - Locator = NA
Column = 14 - Name = Z.Simu.Y.Pb.S6 - Locator = NA
Column = 15 - Name = Z.Simu.Y.Pb.S7 - Locator = NA
Column = 16 - Name = Z.Simu.Y.Pb.S8 - Locator = NA
Column = 17 - Name = Z.Simu.Y.Pb.S9 - Locator = NA
Column = 18 - Name = Z.Simu.Y.Pb.S10 - Locator = NA
Column = 19 - Name = Stats.MINI - Locator = NA
Column = 20 - Name = Stats.MAXI - Locator = NA
Column = 21 - Name = Stats.MEAN - Locator = NA
Column = 22 - Name = Stats.STDV - Locator = z1

We visualize a Conditional Simulation in Raw Scale

In [40]:
if graphics:
    gp.plot(mygrid, "Z.Simu.Y.Pb.S1")
    gp.plot(mydb, nameColor="Pb")
    gp.decoration(title="One simulation of Pb in Raw Scale")
No description has been provided for this image

Let us now average the conditional simulations in order to have a comparison with the estimation by kriging.

In [41]:
mygrid.deleteColumn("Stats.*")
mygrid.statisticsBySample(["Z.Simu.*"], [gl.EStatOption.MEAN], True)
if verbose:
    dbfmt = gl.DbStringFormat()
    dbfmt.setFlags(flag_resume=True, flag_stats=True)
    mygrid.display(dbfmt)
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 20
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Data Base Statistics
--------------------
1 - Name x1 - Locator x1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     107.850
 Maximum value       =     144.850
 Mean value          =     126.350
 Standard Deviation  =      10.824
 Variance            =     117.167
2 - Name x2 - Locator x2
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     481.660
 Maximum value       =     515.160
 Mean value          =     498.410
 Standard Deviation  =       9.814
 Variance            =      96.313
3 - Name Neigh.Pb.Number - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       6.000
 Mean value          =       5.381
 Standard Deviation  =       1.538
 Variance            =       2.366
4 - Name Neigh.Pb.MaxDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.935
 Maximum value       =       9.999
 Mean value          =       5.598
 Standard Deviation  =       2.497
 Variance            =       6.233
5 - Name Neigh.Pb.MinDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.014
 Maximum value       =       9.978
 Mean value          =       3.484
 Standard Deviation  =       2.514
 Variance            =       6.320
6 - Name Neigh.Pb.NbNESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       1.000
 Mean value          =       1.000
 Standard Deviation  =       0.000
 Variance            =       0.000
7 - Name Neigh.Pb.NbCESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.000
 Maximum value       =       0.000
 Mean value          =       0.000
 Standard Deviation  =       0.000
 Variance            =       0.000
8 - Name Kriging.Pb.estim - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.408
 Maximum value       =      11.493
 Mean value          =       6.114
 Standard Deviation  =       0.641
 Variance            =       0.410
9 - Name Kriging.Pb.stdev - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       0.234
 Maximum value       =       1.655
 Mean value          =       1.540
 Standard Deviation  =       0.160
 Variance            =       0.026
10 - Name Z.Simu.Y.Pb.S1 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.593
 Standard Deviation  =       2.989
 Variance            =       8.931
11 - Name Z.Simu.Y.Pb.S2 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.889
 Standard Deviation  =       2.596
 Variance            =       6.741
12 - Name Z.Simu.Y.Pb.S3 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.252
 Standard Deviation  =       2.771
 Variance            =       7.676
13 - Name Z.Simu.Y.Pb.S4 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.777
 Standard Deviation  =       2.993
 Variance            =       8.957
14 - Name Z.Simu.Y.Pb.S5 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.718
 Standard Deviation  =       2.901
 Variance            =       8.415
15 - Name Z.Simu.Y.Pb.S6 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.600
 Standard Deviation  =       2.911
 Variance            =       8.474
16 - Name Z.Simu.Y.Pb.S7 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.665
 Standard Deviation  =       2.893
 Variance            =       8.371
17 - Name Z.Simu.Y.Pb.S8 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.545
 Standard Deviation  =       2.842
 Variance            =       8.079
18 - Name Z.Simu.Y.Pb.S9 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.700
 Standard Deviation  =       2.870
 Variance            =       8.239
19 - Name Z.Simu.Y.Pb.S10 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       6.088
 Standard Deviation  =       2.727
 Variance            =       7.439
20 - Name Stats.MEAN - Locator z1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.340
 Maximum value       =      10.577
 Mean value          =       6.483
 Standard Deviation  =       1.017
 Variance            =       1.035

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = NA
Column = 7 - Name = Kriging.Pb.estim - Locator = NA
Column = 8 - Name = Kriging.Pb.stdev - Locator = NA
Column = 9 - Name = Z.Simu.Y.Pb.S1 - Locator = NA
Column = 10 - Name = Z.Simu.Y.Pb.S2 - Locator = NA
Column = 11 - Name = Z.Simu.Y.Pb.S3 - Locator = NA
Column = 12 - Name = Z.Simu.Y.Pb.S4 - Locator = NA
Column = 13 - Name = Z.Simu.Y.Pb.S5 - Locator = NA
Column = 14 - Name = Z.Simu.Y.Pb.S6 - Locator = NA
Column = 15 - Name = Z.Simu.Y.Pb.S7 - Locator = NA
Column = 16 - Name = Z.Simu.Y.Pb.S8 - Locator = NA
Column = 17 - Name = Z.Simu.Y.Pb.S9 - Locator = NA
Column = 18 - Name = Z.Simu.Y.Pb.S10 - Locator = NA
Column = 19 - Name = Stats.MEAN - Locator = z1

Displaying the average of the Conditional Simulations

In [42]:
if graphics:
    gp.plot(mygrid, "Stats*MEAN")
    gp.plot(mydb, nameColor="Pb")
    gp.decoration(title="Mean of Pb simulations")
No description has been provided for this image

Multivariate case¶

The Gaussian transform of the Pb variable has already been calculated. It suffices to perform the Gaussian transform of the Zn variable

In [43]:
mydb.setLocator("Zn", gl.ELoc.Z)
myanamZn = gl.AnamHermite(30)
myanamZn.fitFromLocator(mydb)
if verbose:
    myanamZn.display()
Hermitian Anamorphosis
----------------------
Minimum absolute value for Y  = -2.5
Maximum absolute value for Y  = 2.6
Minimum absolute value for Z  = 1.1469
Maximum absolute value for Z  = 12.1276
Minimum practical value for Y = -2.5
Maximum practical value for Y = 2.6
Minimum practical value for Z = 1.1469
Maximum practical value for Z = 12.1276
Mean                          = 2.88061
Variance                      = 2.76263
Number of Hermite polynomials = 30
Normalized coefficients for Hermite polynomials (punctual variable)
                 [,  0]     [,  1]     [,  2]     [,  3]     [,  4]     [,  5]     [,  6]
      [  0,]      2.881     -1.277      0.877     -0.447     -0.095      0.294     -0.121
      [  7,]     -0.087      0.134     -0.029     -0.087      0.069      0.034     -0.065
      [ 14,]      0.005      0.044     -0.026     -0.020      0.034      0.001     -0.033
      [ 21,]      0.010      0.027     -0.016     -0.019      0.016      0.012     -0.014
      [ 28,]     -0.005      0.011
In [44]:
if graphics:
    gp.plot(myanamZn)
No description has been provided for this image

We convert the raw data into its Gaussian equivalent

In [45]:
mydb.setLocator("Zn", gl.ELoc.Z)
err = myanamZn.rawToGaussianByLocator(mydb)
if verbose:
    mydb.display()
Data Base Characteristics
=========================

Data Base Summary
-----------------
File is organized as a set of isolated points
Space dimension              = 2
Number of Columns            = 10
Total number of samples      = 102
Number of active samples     = 99

Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = X - Locator = x1
Column = 2 - Name = Y - Locator = x2
Column = 3 - Name = Zn - Locator = NA
Column = 4 - Name = Pb - Locator = NA
Column = 5 - Name = sel - Locator = sel
Column = 6 - Name = Xvalid.Pb.esterr - Locator = NA
Column = 7 - Name = Xvalid.Pb.stderr - Locator = NA
Column = 8 - Name = Y.Pb - Locator = NA
Column = 9 - Name = Y.Zn - Locator = z1

We now perform the multivariate variogram calculation

In [46]:
mydb.setLocators(["Y.Pb", "Y.Zn"], gl.ELoc.Z)
myvario = gl.Vario(myvarioParam)
err = myvario.compute(mydb, gl.ECalcVario.VARIOGRAM)
mymodelM = gl.Model.createFromDb(mydb)
err = mymodelM.fit(myvario, [gl.ECov.EXPONENTIAL])
if graphics:
    axs = gp.varmod(myvario, mymodelM)
    gp.decoration(title="Multivariate Model")
No description has been provided for this image

We perform 10 bivariate conditional simulations (deleting the previous monovariate simulation outcomes first for better legibility)

In [47]:
mygrid.deleteColumn("Z.Simu*")
err = gl.simtub(mydb, mygrid, mymodelM, myneigh, nbsimu=10)
if verbose:
    mygrid.display()
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 30
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = NA
Column = 7 - Name = Kriging.Pb.estim - Locator = NA
Column = 8 - Name = Kriging.Pb.stdev - Locator = NA
Column = 9 - Name = Stats.MEAN - Locator = NA
Column = 10 - Name = Simu.Y.Pb.S1 - Locator = z1
Column = 11 - Name = Simu.Y.Pb.S2 - Locator = z2
Column = 12 - Name = Simu.Y.Pb.S3 - Locator = z3
Column = 13 - Name = Simu.Y.Pb.S4 - Locator = z4
Column = 14 - Name = Simu.Y.Pb.S5 - Locator = z5
Column = 15 - Name = Simu.Y.Pb.S6 - Locator = z6
Column = 16 - Name = Simu.Y.Pb.S7 - Locator = z7
Column = 17 - Name = Simu.Y.Pb.S8 - Locator = z8
Column = 18 - Name = Simu.Y.Pb.S9 - Locator = z9
Column = 19 - Name = Simu.Y.Pb.S10 - Locator = z10
Column = 20 - Name = Simu.Y.Zn.S1 - Locator = z11
Column = 21 - Name = Simu.Y.Zn.S2 - Locator = z12
Column = 22 - Name = Simu.Y.Zn.S3 - Locator = z13
Column = 23 - Name = Simu.Y.Zn.S4 - Locator = z14
Column = 24 - Name = Simu.Y.Zn.S5 - Locator = z15
Column = 25 - Name = Simu.Y.Zn.S6 - Locator = z16
Column = 26 - Name = Simu.Y.Zn.S7 - Locator = z17
Column = 27 - Name = Simu.Y.Zn.S8 - Locator = z18
Column = 28 - Name = Simu.Y.Zn.S9 - Locator = z19
Column = 29 - Name = Simu.Y.Zn.S10 - Locator = z20

We back-transform each set of simulation outcomes using its own Gaussian Anamorphosis function. Finally we delete the Gaussian variables and ask for the statistics on the simulated variables in the Raw Scale.

In [48]:
err = myanamZn.gaussianToRaw(mygrid, "Simu.Y.Zn*")
err = myanamPb.gaussianToRaw(mygrid, "Simu.Y.Pb*")
mygrid.deleteColumn("Simu.Y*")
mygrid.deleteColumn("Stats.*")
mygrid.statisticsBySample(
    ["Z.Simu.*"],
    [
        gl.EStatOption.MINI,
        gl.EStatOption.MAXI,
        gl.EStatOption.MEAN,
        gl.EStatOption.STDV,
    ],
    True,
)
if verbose:
    dbfmt = gl.DbStringFormat()
    dbfmt.setFlags(flag_resume=True, flag_stats=True)
    mygrid.display(dbfmt)
Data Base Grid Characteristics
==============================

Data Base Summary
-----------------
File is organized as a regular grid
Space dimension              = 2
Number of Columns            = 33
Total number of samples      = 5100

Grid characteristics:
---------------------
Origin :     107.850    481.660
Mesh   :       0.500      0.500
Number :          75         68

Data Base Statistics
--------------------
1 - Name x1 - Locator x1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     107.850
 Maximum value       =     144.850
 Mean value          =     126.350
 Standard Deviation  =      10.824
 Variance            =     117.167
2 - Name x2 - Locator x2
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =     481.660
 Maximum value       =     515.160
 Mean value          =     498.410
 Standard Deviation  =       9.814
 Variance            =      96.313
3 - Name Neigh.Pb.Number - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       6.000
 Mean value          =       5.381
 Standard Deviation  =       1.538
 Variance            =       2.366
4 - Name Neigh.Pb.MaxDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.935
 Maximum value       =       9.999
 Mean value          =       5.598
 Standard Deviation  =       2.497
 Variance            =       6.233
5 - Name Neigh.Pb.MinDist - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.014
 Maximum value       =       9.978
 Mean value          =       3.484
 Standard Deviation  =       2.514
 Variance            =       6.320
6 - Name Neigh.Pb.NbNESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       1.000
 Maximum value       =       1.000
 Mean value          =       1.000
 Standard Deviation  =       0.000
 Variance            =       0.000
7 - Name Neigh.Pb.NbCESect - Locator NA
 Nb of data          =        5100
 Nb of active values =        4596
 Minimum value       =       0.000
 Maximum value       =       0.000
 Mean value          =       0.000
 Standard Deviation  =       0.000
 Variance            =       0.000
8 - Name Kriging.Pb.estim - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.408
 Maximum value       =      11.493
 Mean value          =       6.114
 Standard Deviation  =       0.641
 Variance            =       0.410
9 - Name Kriging.Pb.stdev - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       0.234
 Maximum value       =       1.655
 Mean value          =       1.540
 Standard Deviation  =       0.160
 Variance            =       0.026
10 - Name Z.Simu.Y.Zn.S1 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.775
 Standard Deviation  =       1.509
 Variance            =       2.276
11 - Name Z.Simu.Y.Zn.S2 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.931
 Standard Deviation  =       1.737
 Variance            =       3.016
12 - Name Z.Simu.Y.Zn.S3 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.984
 Standard Deviation  =       1.780
 Variance            =       3.169
13 - Name Z.Simu.Y.Zn.S4 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.830
 Standard Deviation  =       1.567
 Variance            =       2.455
14 - Name Z.Simu.Y.Zn.S5 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.845
 Standard Deviation  =       1.598
 Variance            =       2.555
15 - Name Z.Simu.Y.Zn.S6 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.883
 Standard Deviation  =       1.752
 Variance            =       3.071
16 - Name Z.Simu.Y.Zn.S7 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.991
 Standard Deviation  =       1.769
 Variance            =       3.130
17 - Name Z.Simu.Y.Zn.S8 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.871
 Standard Deviation  =       1.660
 Variance            =       2.756
18 - Name Z.Simu.Y.Zn.S9 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.925
 Standard Deviation  =       1.686
 Variance            =       2.843
19 - Name Z.Simu.Y.Zn.S10 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =      12.128
 Mean value          =       2.829
 Standard Deviation  =       1.511
 Variance            =       2.283
20 - Name Z.Simu.Y.Pb.S1 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.724
 Standard Deviation  =       1.762
 Variance            =       3.106
21 - Name Z.Simu.Y.Pb.S2 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.793
 Standard Deviation  =       1.811
 Variance            =       3.281
22 - Name Z.Simu.Y.Pb.S3 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.833
 Standard Deviation  =       1.721
 Variance            =       2.961
23 - Name Z.Simu.Y.Pb.S4 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.842
 Standard Deviation  =       1.732
 Variance            =       2.999
24 - Name Z.Simu.Y.Pb.S5 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.730
 Standard Deviation  =       1.784
 Variance            =       3.183
25 - Name Z.Simu.Y.Pb.S6 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.535
 Standard Deviation  =       1.602
 Variance            =       2.566
26 - Name Z.Simu.Y.Pb.S7 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.788
 Standard Deviation  =       1.783
 Variance            =       3.179
27 - Name Z.Simu.Y.Pb.S8 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.848
 Standard Deviation  =       1.941
 Variance            =       3.767
28 - Name Z.Simu.Y.Pb.S9 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.887
 Standard Deviation  =       1.814
 Variance            =       3.291
29 - Name Z.Simu.Y.Pb.S10 - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.003
 Maximum value       =      12.978
 Mean value          =       5.840
 Standard Deviation  =       1.791
 Variance            =       3.208
30 - Name Stats.MINI - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       1.147
 Maximum value       =       4.178
 Mean value          =       1.835
 Standard Deviation  =       0.286
 Variance            =       0.082
31 - Name Stats.MAXI - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       3.445
 Maximum value       =      12.978
 Mean value          =       9.193
 Standard Deviation  =       2.062
 Variance            =       4.252
32 - Name Stats.MEAN - Locator NA
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       2.376
 Maximum value       =       9.400
 Mean value          =       4.334
 Standard Deviation  =       0.595
 Variance            =       0.354
33 - Name Stats.STDV - Locator z1
 Nb of data          =        5100
 Nb of active values =        5100
 Minimum value       =       0.579
 Maximum value       =       3.887
 Mean value          =       2.116
 Standard Deviation  =       0.477
 Variance            =       0.227

Variables
---------
Column = 0 - Name = x1 - Locator = x1
Column = 1 - Name = x2 - Locator = x2
Column = 2 - Name = Neigh.Pb.Number - Locator = NA
Column = 3 - Name = Neigh.Pb.MaxDist - Locator = NA
Column = 4 - Name = Neigh.Pb.MinDist - Locator = NA
Column = 5 - Name = Neigh.Pb.NbNESect - Locator = NA
Column = 6 - Name = Neigh.Pb.NbCESect - Locator = NA
Column = 7 - Name = Kriging.Pb.estim - Locator = NA
Column = 8 - Name = Kriging.Pb.stdev - Locator = NA
Column = 9 - Name = Z.Simu.Y.Zn.S1 - Locator = NA
Column = 10 - Name = Z.Simu.Y.Zn.S2 - Locator = NA
Column = 11 - Name = Z.Simu.Y.Zn.S3 - Locator = NA
Column = 12 - Name = Z.Simu.Y.Zn.S4 - Locator = NA
Column = 13 - Name = Z.Simu.Y.Zn.S5 - Locator = NA
Column = 14 - Name = Z.Simu.Y.Zn.S6 - Locator = NA
Column = 15 - Name = Z.Simu.Y.Zn.S7 - Locator = NA
Column = 16 - Name = Z.Simu.Y.Zn.S8 - Locator = NA
Column = 17 - Name = Z.Simu.Y.Zn.S9 - Locator = NA
Column = 18 - Name = Z.Simu.Y.Zn.S10 - Locator = NA
Column = 19 - Name = Z.Simu.Y.Pb.S1 - Locator = NA
Column = 20 - Name = Z.Simu.Y.Pb.S2 - Locator = NA
Column = 21 - Name = Z.Simu.Y.Pb.S3 - Locator = NA
Column = 22 - Name = Z.Simu.Y.Pb.S4 - Locator = NA
Column = 23 - Name = Z.Simu.Y.Pb.S5 - Locator = NA
Column = 24 - Name = Z.Simu.Y.Pb.S6 - Locator = NA
Column = 25 - Name = Z.Simu.Y.Pb.S7 - Locator = NA
Column = 26 - Name = Z.Simu.Y.Pb.S8 - Locator = NA
Column = 27 - Name = Z.Simu.Y.Pb.S9 - Locator = NA
Column = 28 - Name = Z.Simu.Y.Pb.S10 - Locator = NA
Column = 29 - Name = Stats.MINI - Locator = NA
Column = 30 - Name = Stats.MAXI - Locator = NA
Column = 31 - Name = Stats.MEAN - Locator = NA
Column = 32 - Name = Stats.STDV - Locator = z1

Categorical Variable¶

We compare the initial variable 'Pb' with a set of disjoint intervals. The 'Pb' values varying from 3 to 12.7, we consider three classes:

  • values below 4
  • values between 4 and 6
  • values above 6

We first build the indicators for each class

In [49]:
limits = gl.Limits([np.nan, 4.0, 6.0, np.nan])
if verbose:
    limits.display()
Bound( 1 ) : ] -Inf ; 4 [
Bound( 2 ) : [ 4 ; 6 [
Bound( 3 ) : [ 6 ;  +Inf [

We apply the set of limits previously defined in order to transform the input variable into Indicators of the different classes.

In [50]:
err = limits.toIndicator(mydb, "Pb")
if verbose:
    mydb.display()
Data Base Characteristics
=========================

Data Base Summary
-----------------
File is organized as a set of isolated points
Space dimension              = 2
Number of Columns            = 13
Total number of samples      = 102
Number of active samples     = 99

Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = X - Locator = x1
Column = 2 - Name = Y - Locator = x2
Column = 3 - Name = Zn - Locator = NA
Column = 4 - Name = Pb - Locator = NA
Column = 5 - Name = sel - Locator = sel
Column = 6 - Name = Xvalid.Pb.esterr - Locator = NA
Column = 7 - Name = Xvalid.Pb.stderr - Locator = NA
Column = 8 - Name = Y.Pb - Locator = NA
Column = 9 - Name = Y.Zn - Locator = NA
Column = 10 - Name = Indicator.Pb.Class.S1 - Locator = z1
Column = 11 - Name = Indicator.Pb.Class.S2 - Locator = z2
Column = 12 - Name = Indicator.Pb.Class.S3 - Locator = z3

We calculate the variogram of the Indicators for future use

In [51]:
myvarioindParam = gl.VarioParam()
myvarioindParam.addDir(mydir)
myvarioInd = gl.Vario(myvarioindParam)
err = myvarioInd.compute(mydb, gl.ECalcVario.VARIOGRAM)
if verbose:
    myvarioInd.display()
Variogram characteristics
=========================
Number of variable(s)       = 3
Number of direction(s)      = 1
Space dimension             = 2
Variable(s)                 = [Indicator.Pb.Class.S1 Indicator.Pb.Class.S2 Indicator.Pb.Class.S3]

Variance-Covariance Matrix
                 [,  0]     [,  1]     [,  2]
      [  0,]      0.107     -0.062     -0.044
      [  1,]     -0.062      0.250     -0.187
      [  2,]     -0.044     -0.187      0.231

Direction #1
------------
Number of lags              = 10
Direction coefficients      =       1.000      0.000
Direction angles (degrees)  =       0.000
Tolerance on direction      =      90.000 (degrees)
Calculation lag             =       1.000
Tolerance on distance       =      50.000 (Percent of the lag value)

For variable 1
       Rank     Npairs   Distance      Value
          0      3.000      0.389      0.000
          1    123.000      1.081      0.081
          2    183.000      2.038      0.126
          3    205.000      3.006      0.156
          4    231.000      4.013      0.132
          5    229.000      5.036      0.159
          6    198.000      5.962      0.152
          7    187.000      7.000      0.107
          8    204.000      7.996      0.096
          9    184.000      8.990      0.068

For variables 2 and 1
       Rank     Npairs   Distance      Value
          0      3.000      0.389      0.000
          1    123.000      1.081     -0.065
          2    183.000      2.038     -0.077
          3    205.000      3.006     -0.085
          4    231.000      4.013     -0.093
          5    229.000      5.036     -0.085
          6    198.000      5.962     -0.061
          7    187.000      7.000     -0.045
          8    204.000      7.996     -0.042
          9    184.000      8.990     -0.038

For variable 2
       Rank     Npairs   Distance      Value
          0      3.000      0.389      0.167
          1    123.000      1.081      0.199
          2    183.000      2.038      0.221
          3    205.000      3.006      0.251
          4    231.000      4.013      0.292
          5    229.000      5.036      0.258
          6    198.000      5.962      0.237
          7    187.000      7.000      0.254
          8    204.000      7.996      0.228
          9    184.000      8.990      0.234

For variables 3 and 1
       Rank     Npairs   Distance      Value
          0      3.000      0.389      0.000
          1    123.000      1.081     -0.016
          2    183.000      2.038     -0.049
          3    205.000      3.006     -0.071
          4    231.000      4.013     -0.039
          5    229.000      5.036     -0.074
          6    198.000      5.962     -0.091
          7    187.000      7.000     -0.061
          8    204.000      7.996     -0.054
          9    184.000      8.990     -0.030

For variables 3 and 2
       Rank     Npairs   Distance      Value
          0      3.000      0.389     -0.167
          1    123.000      1.081     -0.134
          2    183.000      2.038     -0.145
          3    205.000      3.006     -0.166
          4    231.000      4.013     -0.199
          5    229.000      5.036     -0.172
          6    198.000      5.962     -0.177
          7    187.000      7.000     -0.209
          8    204.000      7.996     -0.186
          9    184.000      8.990     -0.196

For variable 3
       Rank     Npairs   Distance      Value
          0      3.000      0.389      0.167
          1    123.000      1.081      0.150
          2    183.000      2.038      0.194
          3    205.000      3.006      0.237
          4    231.000      4.013      0.238
          5    229.000      5.036      0.247
          6    198.000      5.962      0.268
          7    187.000      7.000      0.270
          8    204.000      7.996      0.240
          9    184.000      8.990      0.226
In [52]:
ax = gp.varmod(myvarioInd)
No description has been provided for this image

Then we build a categorical variable which gives the index of the class to which each sample belongs

In [53]:
err = limits.toCategory(mydb, "Pb")
if verbose:
    dbfmt = gl.DbStringFormat()
    dbfmt.setFlags(flag_stats=True)
    dbfmt.setNames(["Category*"])
    dbfmt.setMode(2)
    mydb.display(dbfmt)
Data Base Characteristics
=========================

Data Base Summary
-----------------
File is organized as a set of isolated points
Space dimension              = 2
Number of Columns            = 14
Total number of samples      = 102
Number of active samples     = 99

Data Base Statistics
--------------------
14 - Name Category.Pb - Locator z1
 Nb of data          =         102
 Nb of active values =          99
 Class          1 =          12 (     12.121%)
 Class          2 =          51 (     51.515%)
 Class          3 =          36 (     36.364%)

Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = X - Locator = x1
Column = 2 - Name = Y - Locator = x2
Column = 3 - Name = Zn - Locator = NA
Column = 4 - Name = Pb - Locator = NA
Column = 5 - Name = sel - Locator = sel
Column = 6 - Name = Xvalid.Pb.esterr - Locator = NA
Column = 7 - Name = Xvalid.Pb.stderr - Locator = NA
Column = 8 - Name = Y.Pb - Locator = NA
Column = 9 - Name = Y.Zn - Locator = NA
Column = 10 - Name = Indicator.Pb.Class.S1 - Locator = NA
Column = 11 - Name = Indicator.Pb.Class.S2 - Locator = NA
Column = 12 - Name = Indicator.Pb.Class.S3 - Locator = NA
Column = 13 - Name = Category.Pb - Locator = z1
In [ ]: