variography¶
Variography¶
%matplotlib inline
# loading packages
import pandas as pd
import numpy as np
geolime package can be imported as such:
import geolime as glm
Data¶
A dataset is available for testing. Here we load some geological data and display the first five rows.
assay_2d_tnb = "../../tests/data/variography/assay_2d_tnb.csv"
data = pd.read_csv(assay_2d_tnb)
data.head()
| Unnamed: 0 | holeid | x | y | z | from | to | ni | co | fe | mgo | sio2 | domain | layer | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0 | 03PD2018L370O150 ... | 240242.692 | 383322.588 | 282.916 | 0.0 | 1.0 | 2.768 | 0.019 | 7.84 | 26.88 | 44.76 | D | SAP |
| 1 | 1 | 03PD2018L370O270 ... | 240250.718 | 383440.213 | 286.957 | 0.0 | 1.0 | 0.699 | 0.021 | 9.73 | 28.16 | 41.88 | D | SAP |
| 2 | 2 | 03PD2018L404O430 ... | 240280.549 | 383602.858 | 289.022 | 0.0 | 1.0 | 0.440 | 0.020 | 8.12 | 28.93 | 43.47 | D | SAP |
| 3 | 3 | 03PD2018L410O310 ... | 240291.812 | 383480.851 | 305.434 | 0.0 | 1.0 | 1.418 | 0.027 | 7.02 | 29.75 | 47.01 | D | SAP |
| 4 | 4 | 03PD2018L410O350 ... | 240288.891 | 383520.723 | 302.360 | 1.2 | 2.0 | 1.560 | 0.026 | 11.84 | 20.40 | 47.19 | T | SAP |
Experimental semivariogram¶
Prepare arrays to be used as input for experimental semivariogram computation.
# extract coordinates and values
data = pd.read_csv(assay_2d_tnb)
coords_df = data[['x', 'y']]
coords = coords_df.to_numpy()
values_df = data[['ni']]
values = values_df.to_numpy()
Define input lags and compute the experimental semivariogram.
# define input lags
input_lags = [0, 13.34, 26.68, 40.02, 53.36, 66.7,
80.04, 93.38, 106.72, 120.06]
# compute variogram
vario_exp_az_north = glm.variogram(coords=coords, values=values, lags=input_lags,
tol=6.67, geographic_azimuth=0, dip=0, pitch=0,
atol=45, slice_width=np.nan, slice_height=np.nan,
method="semivariogram")
# compute variogram
vario_exp_az_east = glm.variogram(coords=coords, values=values, lags=input_lags,
tol=6.67, geographic_azimuth=90, dip=0, pitch=0,
atol=45, slice_width=np.nan, slice_height=np.nan,
method="semivariogram")
The resulting dataframe contains the following columns
vario_exp_az_north.head()
| rank | npairs | lag | avgdist | vario | indices | |
|---|---|---|---|---|---|---|
| 0 | 0 | 12 | 0.00 | 5.97538 | 0.537738 | [[17, 51], [24, 54], [32, 104], [34, 118], [35... |
| 1 | 1 | 597 | 13.34 | 13.2226 | 0.506487 | [[0, 259], [1, 83], [1, 92], [2, 269], [7, 268... |
| 2 | 2 | 916 | 26.68 | 26.6512 | 0.555162 | [[0, 255], [0, 260], [1, 76], [1, 93], [1, 98]... |
| 3 | 3 | 999 | 40.02 | 40.2595 | 0.610749 | [[0, 29], [0, 104], [1, 89], [1, 90], [1, 99],... |
| 4 | 4 | 955 | 53.36 | 53.2865 | 0.716520 | [[0, 32], [0, 91], [0, 95], [0, 242], [1, 25],... |
vario_exp_az_east.head()
| rank | npairs | lag | avgdist | vario | indices | |
|---|---|---|---|---|---|---|
| 0 | 0 | 1 | 0.00 | 6.23131 | 0.004050 | [[203, 207]] |
| 1 | 1 | 603 | 13.34 | 13.4262 | 0.475794 | [[0, 254], [2, 46], [5, 266], [9, 263], [10, 2... |
| 2 | 2 | 854 | 26.68 | 26.9149 | 0.595443 | [[0, 249], [1, 26], [1, 96], [1, 97], [1, 251]... |
| 3 | 3 | 1072 | 40.02 | 40.2483 | 0.647403 | [[0, 112], [0, 241], [1, 105], [1, 244], [1, 2... |
| 4 | 4 | 1107 | 53.36 | 53.338 | 0.691535 | [[0, 33], [0, 111], [0, 119], [0, 235], [1, 3]... |
vario_exp = [vario_exp_az_north, vario_exp_az_east]
It also contains useful metadata
vario_exp_az_north.attrs
{'method': 'semivariogram',
'geographic_azimuth': 0,
'dip': 0,
'pitch': 0,
'tol': 6.67,
'atol': 45.0,
'trend': 90.0,
'plunge': 0.0,
'variable_mean': 1.8995185185185182,
'variable_var': 0.7367737459533608,
'number_of_lags': 10,
'input_lags': [0,
13.34,
26.68,
40.02,
53.36,
66.7,
80.04,
93.38,
106.72,
120.06],
'slice_width': nan,
'slice_height': nan,
'dim': 2,
'directional_vector': array([6.123234e-17, 1.000000e+00])}
We can plot the experimental semivariogram dataframe with the following command:
glm.plot_semivariogram(variograms=vario_exp, model=None, display_sill=True)
(<Figure size 432x288 with 1 Axes>, <AxesSubplot:>)
Variogram model¶
Experimental semivariogram can be fitted the following way. We first need to define a target grid.
# create a grid
target_grid = glm.create_grid_2d(origin=[0, 0], n_dim=[30, 30], cell_size=[5, 5], angle=[0.0, 0.0, 0.0])
We also need to define a covariance model for fitting
# define a covariance model
cov = glm.Nugget() + 0.5 * glm.Spherical(2)
cov.print()
Model with 2 components -------------------------------------------------- Component 1 : Variance : 1 Covariance type : Nugget --------- Component 2 : Variance : 0.5 Covariance type : Spherical Scales : [1. 1.] Angles : 0.0 Rotation matrix [[ 1. -0.] [ 0. 1.]] --------- Total sill = 1.5
We can fit the experimental variogram
# fit model
cov_model = glm.model_fit(variograms=vario_exp, cov=cov, constraint=None)
glm.plot_semivariogram(variograms=[vario_exp_az_north], model=cov_model)
(<Figure size 432x288 with 1 Axes>, <AxesSubplot:>)
cov_model.print()
Model with 2 components -------------------------------------------------- Component 1 : Variance : 0.40751929468954295 Covariance type : Nugget --------- Component 2 : Variance : 0.3432202933243947 Covariance type : Spherical Scales : [79.52391655 83.24973102] Angles : 14.651468428093645 Rotation matrix [[ 0.96748235 -0.25293854] [ 0.25293854 0.96748235]] --------- Total sill = 0.7507395880139376
We apply model regularization with the following call:
pattern_generator = glm.GridPatternGenerator(target_grid, [5, 5])
mr = glm.model_regularize(vario_exp_az_north, cov, pattern_generator, nlag=10, border_ratio=1, normalize=True)
mr
| lag | vario | npairs | |
|---|---|---|---|
| 0 | 0.00 | 0.000000 | 1.0 |
| 1 | 13.34 | 0.203300 | 1.0 |
| 2 | 26.68 | 0.440564 | 1.0 |
| 3 | 40.02 | 0.651105 | 1.0 |
| 4 | 53.36 | 0.823227 | 1.0 |
| 5 | 66.70 | 0.943547 | 1.0 |
| 6 | 80.04 | 0.997068 | 1.0 |
| 7 | 93.38 | 1.000000 | 1.0 |
| 8 | 106.72 | 1.000000 | 1.0 |
| 9 | 120.06 | 1.000000 | 1.0 |
Finally, the experimental variogram can be fitted and plotted with the following call:
# fit regularized model
rmodel = glm.model_fit(variograms=[mr], cov=cov, constraint=None)
glm.plot_semivariogram(variograms=[mr], model=rmodel)
(<Figure size 432x288 with 1 Axes>, <AxesSubplot:>)