🎯 But pédagogique¶
Comprendre l’influence du choix du modèle de variogramme (ou covariance) sur l’estimation de la variance des blocs et la dispersion des valeurs. L’objectif est de montrer que, pour une même variance globale, différents modèles peuvent engendrer des structures de variabilité spatiale très différentes, ce qui affecte directement :
- l’estimation des ressources sur des blocs de taille donnée,
- la fiabilité des estimations dans les domaines peu ou mal échantillonnés,
- le comportement de l’estimateur (krigeage, IDW, etc.) face à l’échelle du problème.
🔬 Activité interactive¶
À l’aide des sliders, vous pouvez :
- Modifier les paramètres des modèles (nugget, sill, portée).
- Comparer plusieurs modèles sur :
- le graphique des variogrammes,
- la courbe de variance de bloc,
- la dispersion (différence entre blocs petits et blocs de référence).
- Analyser les paramètres du variogramme peut produire des résultats très différents sur la dispersion spatiale.
❓ Questions de réflexion¶
- Pourquoi deux modèles avec le même plateau donnent-ils des variances de bloc différentes ?
- Que traduit une dispersion élevée à petite échelle ?
- Quel modèle serait le plus pénalisant pour l’estimation des ressources dans un gisement mal foré ?
- Que se passe-t-il si on sous-estime ou surestime la portée principale du variogramme ?
Source
import numpy as np
import matplotlib.pyplot as plt
from numpy.polynomial.legendre import leggauss
import ipywidgets as widgets
from IPython.display import display
# --- Style graphique global ---
plt.rcParams.update({
'font.family': 'DejaVu Sans',
'font.size': 13,
'axes.titlesize': 15,
'axes.labelsize': 13,
'xtick.labelsize': 11,
'ytick.labelsize': 11,
'legend.fontsize': 11,
'lines.linewidth': 2,
'axes.grid': True,
'grid.linestyle': '--',
'grid.alpha': 0.6,
'figure.facecolor': 'white',
'axes.facecolor': '#f9f9f9'
})
# --- Covariance sphérique ---
def spherical_covariance(h, range_, sill):
cov = np.zeros_like(h)
mask = h <= range_
hr = h[mask] / range_
cov[mask] = sill * (1 - 1.5 * hr + 0.5 * hr**3)
return cov
def spherical_variogram(h, nugget, sill, range_):
return nugget + (sill - nugget) * (1.5 * (h / range_) - 0.5 * (h / range_)**3) * (h <= range_) + (sill - nugget) * (h > range_)
def nested_spherical_variogram(h, nugget, sill1, range1, sill2, range2):
gamma1 = spherical_variogram(h, 0, sill1, range1)
gamma2 = spherical_variogram(h, 0, sill2, range2)
return nugget + gamma1 + gamma2
# --- Covariance imbriquée pour la variance de bloc ---
def block_variance_surface(lx, ly, model, n_points=5):
pts_1D, w_1D = leggauss(n_points)
pts_1D = 0.5 * (pts_1D + 1)
w_1D = 0.5 * w_1D
X, Y = np.meshgrid(pts_1D * lx, pts_1D * ly)
points = np.column_stack([X.ravel(), Y.ravel()])
weights = np.outer(w_1D, w_1D).ravel()
diff = (points[:, None, :] - points[None, :, :])
if model['type'] == 'spherical':
h = np.linalg.norm(diff / model['range'], axis=2)
cov_mat = spherical_covariance(h, 1.0, model['sill'] - model['nugget'])
elif model['type'] == 'nested':
h1 = np.linalg.norm(diff / model['range1'], axis=2)
h2 = np.linalg.norm(diff / model['range2'], axis=2)
cov1 = spherical_covariance(h1, 1.0, model['sill1'])
cov2 = spherical_covariance(h2, 1.0, model['sill2'])
cov_mat = cov1 + cov2
else:
raise NotImplementedError("Modèle non supporté.")
return np.sum(weights[:, None] * weights[None, :] * cov_mat)
# --- Fonction principale avec tracé ---
def plot_variogrammes(c0=0.2, c=0.8, range1=40, l_bloc=20,
c0_3=0.2, sill_3=1.7, range_3=60):
variogrammes = [
{"type": "spherical", "nugget": c0, "sill": c0 + c, "range": range1},
{"type": "nested", "nugget": 0.0, "sill1": c0, "range1": 5, "sill2": c, "range2": range1},
{"type": "spherical", "nugget": c0_3, "sill": sill_3, "range": range_3},
]
h_vals = np.linspace(0, 100, 300)
fig, axs = plt.subplots(1, 3, figsize=(18, 5), constrained_layout=True)
# --- 1. Variogrammes ---
for i, v in enumerate(variogrammes):
if v['type'] == 'spherical':
gamma = spherical_variogram(h_vals, v['nugget'], v['sill'], v['range'])
elif v['type'] == 'nested':
gamma = nested_spherical_variogram(h_vals, v['nugget'], v['sill1'], v['range1'], v['sill2'], v['range2'])
axs[0].plot(h_vals, gamma, label=f"Vario {i+1}")
axs[0].axvline(x=l_bloc, color='red', linestyle='--', linewidth=2, label=f'Bloc {l_bloc}×{l_bloc}')
axs[0].set_title("Variogrammes")
axs[0].set_xlabel("Distance h")
axs[0].set_ylabel("γ(h)")
axs[0].legend(frameon=False)
# --- 2. Variance de bloc ---
l_vals = np.arange(1, 41)
for i, v in enumerate(variogrammes):
var_blocs = [block_variance_surface(l, l, model=v, n_points=5) for l in l_vals]
axs[1].plot(l_vals, var_blocs, label=f"Vario {i+1}")
axs[1].set_title("Variance de bloc vs taille")
axs[1].set_xlabel("Taille du bloc (l)")
axs[1].set_ylabel("Variance")
axs[1].legend(frameon=False)
# --- 3. Dispersion = Var(bloc réf) - Var(bloc l) ---
l_vals_disp = np.arange(1, l_bloc + 1)
for i, v in enumerate(variogrammes):
var_bloc_ref = block_variance_surface(l_bloc, l_bloc, model=v, n_points=5)
var_blocs = [block_variance_surface(l, l, model=v, n_points=5) for l in l_vals_disp]
dispersion = np.array(var_blocs) - var_bloc_ref
axs[2].plot(l_vals_disp, dispersion, label=f"Vario {i+1}")
axs[2].set_title(f"Dispersion : Var(l×l) - Var({l_bloc}×{l_bloc})")
axs[2].set_xlabel("Taille du bloc l")
axs[2].set_ylabel("Dispersion")
axs[2].legend(frameon=False)
plt.show()
# --- Sliders interactifs ---
slider_c0 = widgets.FloatSlider(value=0.2, min=0.0, max=0.5, step=0.05, description='c₀ (v1,v2)')
slider_c = widgets.FloatSlider(value=1.0, min=0.1, max=1.5, step=0.1, description='c (v1,v2)')
slider_range1 = widgets.IntSlider(value=40, min=10, max=100, step=5, description='Portée (v1,v2)')
slider_c0_3 = widgets.FloatSlider(value=0.2, min=0.0, max=0.5, step=0.05, description='c₀ vario3')
slider_sill_3 = widgets.FloatSlider(value=1.7, min=0.5, max=2.0, step=0.1, description='Sill vario3')
slider_range_3 = widgets.IntSlider(value=60, min=10, max=100, step=5, description='Portée vario3')
slider_l_bloc = widgets.IntSlider(value=20, min=5, max=40, step=1, description='Taille bloc réf')
group1 = widgets.VBox([
widgets.HTML("<h4 style='margin-bottom:0'>🔹 Paramètres des variogrammes 1 & 2</h4>"),
slider_c0, slider_c, slider_range1
])
group2 = widgets.VBox([
widgets.HTML("<h4 style='margin-bottom:0'>🔹 Paramètres du variogramme 3</h4>"),
slider_c0_3, slider_sill_3, slider_range_3
])
group3 = widgets.VBox([
widgets.HTML("<h4 style='margin-bottom:0'>🔹 Paramètre du bloc de référence</h4>"),
slider_l_bloc
])
ui = widgets.VBox([group1, group2, group3])
# Lien avec la fonction
out = widgets.interactive_output(
plot_variogrammes,
{
'c0': slider_c0,
'c': slider_c,
'range1': slider_range1,
'c0_3': slider_c0_3,
'sill_3': slider_sill_3,
'range_3': slider_range_3,
'l_bloc': slider_l_bloc
}
)
display(ui, out)
Loading...