Atelier interactif pour comprendre comment les paramètres () et () sont calculés
La valeur de correspond à la proportion, représentée par la zone hachurée dans le graphique de la fonction de densité () ; elle s’obtient aussi comme à partir de la fonction de répartition (non présentée).
La valeur de correspond à la moyenne des teneurs situées dans cette zone hachurée (c’est-à-dire, ) sur le graphique de la fonction de densité.
Note : Les paramètres et dépendent de la distribution statistique des teneurs dans le gisement.
Ici, on suppose une loi lognormale de moyenne et de variance .
Toutefois, toute autre distribution pourrait être utilisée selon les caractéristiques réelles du gisement.
Cela souligne l’importance de disposer d’une quantité suffisante d’informations pour bien caractériser la distribution des teneurs.
Source
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm, lognorm
import ipywidgets as widgets
from IPython.display import display
plt.rcParams.update({'font.size': 14}) # Agrandit globalement les polices
def reserve(moy, s2, c_array):
c_array = np.asarray(c_array)
log_term = np.log(s2 / moy**2 + 1)
b = np.sqrt(log_term)
u = np.log(moy) - 0.5 * log_term
log_mc = np.log(moy / c_array)
n1 = log_mc / b - b / 2
n2 = n1 + b
tc = norm.cdf(n1)
qc = moy * norm.cdf(n2)
mc = np.where(tc > 0, qc / tc, 0)
return tc, qc, mc
def plot_interactive(moy, s2, c_cut):
c_range = np.linspace(0.001, 8, 500)
tc_vals, qc_vals, mc_vals = reserve(moy, s2, c_range)
tc_cut, qc_cut, mc_cut = reserve(moy, s2, np.array([c_cut]))
log_term = np.log(s2 / moy**2 + 1)
sigma = np.sqrt(log_term)
mu = np.log(moy) - 0.5 * log_term
x = np.linspace(0.001, 8, 500)
cdf = lognorm.cdf(x, s=sigma, scale=np.exp(mu))
pdf = lognorm.pdf(x, s=sigma, scale=np.exp(mu))
cdf_cut = lognorm.cdf(c_cut, s=sigma, scale=np.exp(mu))
pdf_cut = lognorm.pdf(c_cut, s=sigma, scale=np.exp(mu))
fig, axs = plt.subplots(2, 2, figsize=(12, 9))
axs = axs.ravel()
axs[0].plot(c_range, tc_vals, label=r'$x_c = P(X > c)$')
axs[0].axvline(c_cut, color='red', linestyle='--')
axs[0].scatter([c_cut], tc_cut, color='red')
axs[0].text(c_cut+0.1, tc_cut + 0.03, f"{tc_cut[0]:.2f}", color='red')
axs[0].set(xlim=(0, 8), ylim=(0, 1), xlabel='Teneur de coupure ($c$)', ylabel='$x_c$',
title=r'Proportion au-dessus de $c$')
axs[0].legend()
axs[0].grid()
axs[1].plot(c_range, mc_vals, label=r'$g_c = E[X|X>c]$')
axs[1].axvline(c_cut, color='red', linestyle='--')
axs[1].scatter([c_cut], mc_cut, color='red')
axs[1].text(c_cut+0.1, mc_cut + 2, f"{mc_cut[0]:.2f}", color='red')
axs[1].set(xlim=(0, 8), ylim=(0, 15), xlabel='Teneur de coupure ($c$)', ylabel='$g_c$',
title='Teneur moyenne conditionnelle')
axs[1].legend()
axs[1].grid()
axs[2].plot(c_range, tc_vals * mc_vals, label=r'$x_c \cdot g_c$')
axs[2].axvline(c_cut, color='red', linestyle='--')
axs[2].scatter([c_cut], [tc_cut[0] * mc_cut[0]], color='red')
axs[2].text(c_cut+0.1, tc_cut[0] * mc_cut[0] + 0.03, f"{tc_cut[0] * mc_cut[0]:.2f}", color='red')
axs[2].set(xlim=(0, 8),
ylim=(0, np.max(tc_vals * mc_vals) * 1.1),
xlabel='Teneur de coupure ($c$)',
ylabel=r'$x_c \cdot g_c$',
title='Teneur moyenne récupérée')
axs[2].legend()
axs[2].grid()
axs[3].plot(x, pdf, label='PDF')
axs[3].fill_between(x, 0, pdf, where=(x > c_cut), color='grey', alpha=0.3, hatch='///')
axs[3].axvline(c_cut, color='red', linestyle='--')
axs[3].scatter([c_cut], [pdf_cut], color='red')
axs[3].text(c_cut+0.1, pdf_cut + 0.02, f"{pdf_cut:.2f}", color='red')
axs[3].set(xlim=(0, 8), ylim=(0, np.max(pdf) + 0.1), xlabel='Teneur de coupure ($c$)', ylabel='Densité',
title='Fonction de densité lognormale')
axs[3].legend()
axs[3].grid()
plt.tight_layout()
plt.show()
# Widgets
w_moy = widgets.FloatSlider(value=1.0, min=0.5, max=2.0, step=0.05, description='Moyenne')
w_s2 = widgets.FloatSlider(value=4.0, min=1.0, max=8.0, step=0.25, description='Variance')
w_c = widgets.FloatSlider(value=0.5, min=0.001, max=8.0, step=0.2, description='Coupure')
ui = widgets.VBox([w_moy, w_s2, w_c])
out = widgets.interactive_output(plot_interactive, {'moy': w_moy, 's2': w_s2, 'c_cut': w_c})
display(ui, out)