Aggiunta dati e versione finale di focale lente
This commit is contained in:
+110
-16
@@ -1,6 +1,6 @@
|
||||
import marimo
|
||||
|
||||
__generated_with = "0.25.0"
|
||||
__generated_with = "0.25.2"
|
||||
app = marimo.App()
|
||||
|
||||
|
||||
@@ -163,7 +163,7 @@ def _(np, pd, stats):
|
||||
return campione, rimossi
|
||||
|
||||
|
||||
def stat_chauv(self, prefissi, err_strumentale):
|
||||
def stat_chauv(self, prefissi, err_strumentale, err_della_media = True):
|
||||
campione_df = self.campione # riferimento al DataFrame originale (non sovrascritto)
|
||||
|
||||
if isinstance(prefissi, str):
|
||||
@@ -185,15 +185,16 @@ def _(np, pd, stats):
|
||||
valori_puliti, rimossi = rimuovi_outlier(valori, err_strum)
|
||||
|
||||
n = len(valori_puliti)
|
||||
s = valori_puliti.std(ddof=1) / np.sqrt(n)
|
||||
s = valori_puliti.std(ddof=1)
|
||||
if err_della_media: s /= np.sqrt(n)
|
||||
u = np.sqrt(s**2 + err_strum**2)
|
||||
|
||||
self.analisi_stat.at[i, prefix] = valori_puliti.mean()
|
||||
self.analisi_stat.at[i, f"u{prefix}_strum"] = err_strum
|
||||
self.analisi_stat.at[i, f"u{prefix}_stat"] = s
|
||||
self.analisi_stat.at[i, f"u{prefix}"] = u
|
||||
self.analisi_stat.at[i, f"n{prefix}"] = n
|
||||
self.analisi_stat.at[i, f"out{prefix}"] = rimossi
|
||||
self.analisi_stat.loc[i, prefix] = valori_puliti.mean()
|
||||
self.analisi_stat.loc[i, f"u{prefix}_strum"] = err_strum
|
||||
self.analisi_stat.loc[i, f"u{prefix}_stat"] = s
|
||||
self.analisi_stat.loc[i, f"u{prefix}"] = u
|
||||
self.analisi_stat.loc[i, f"n{prefix}"] = n
|
||||
self.analisi_stat.loc[i, f"out{prefix}"] = rimossi
|
||||
|
||||
return self
|
||||
|
||||
@@ -201,6 +202,95 @@ def _(np, pd, stats):
|
||||
return
|
||||
|
||||
|
||||
@app.cell(hide_code=True)
|
||||
def _(mo):
|
||||
mo.md(r"""
|
||||
# Media pesata
|
||||
|
||||
input:
|
||||
- array x: valori
|
||||
- array ux: uvalori
|
||||
output:
|
||||
- xm: x medio
|
||||
- uxm: sigma x medio
|
||||
- chi2: $\chi^2$
|
||||
- dof: gradi di libertà
|
||||
- Pchi2: $P(\chi^2)$
|
||||
|
||||
data:
|
||||
- prefix: prefisso del nome della variabile
|
||||
""")
|
||||
return
|
||||
|
||||
|
||||
@app.cell
|
||||
def _(np, stats):
|
||||
def media_pesata(x, ux):
|
||||
x, ux = np.asarray(x, float), np.asarray(ux, float)
|
||||
w = 1 / ux**2
|
||||
xm = np.sum(w * x) / np.sum(w)
|
||||
uxm = 1 / np.sqrt(np.sum(w))
|
||||
chi2 = np.sum(w * (x - xm) ** 2)
|
||||
dof = len(x) - 1
|
||||
return xm, uxm, chi2/dof, dof, stats.chi2.sf(chi2, dof)
|
||||
|
||||
def media_pesata_data(self, prefix):
|
||||
return media_pesata(self.analisi_stat[prefix], self.analisi_stat[f"u{prefix}"])
|
||||
|
||||
Data.media_pesata = media_pesata_data
|
||||
return
|
||||
|
||||
|
||||
@app.cell
|
||||
def _(np, stats):
|
||||
def media_pesata_dersimonian_laird(x, ux):
|
||||
x, ux = np.asarray(x, float), np.asarray(ux, float)
|
||||
n = len(x)
|
||||
dof = n - 1
|
||||
|
||||
if n <= 1:
|
||||
raise ValueError("Il dataset deve contenere almeno due misure.")
|
||||
if np.any(ux <= 0):
|
||||
raise ValueError("Le incertezze (ux) devono essere strettamente maggiori di zero.")
|
||||
|
||||
# --- STEP 1: Calcolo standard (Effetti Fissi) per stimare il Chi2 originale ---
|
||||
w_fixed = 1 / ux**2
|
||||
sum_w_fixed = np.sum(w_fixed)
|
||||
xm_fixed = np.sum(w_fixed * x) / sum_w_fixed
|
||||
chi2_fixed = np.sum(w_fixed * (x - xm_fixed) ** 2)
|
||||
|
||||
# --- STEP 2: Stima della varianza extra (Tau^2) ---
|
||||
# Formula di DerSimonian-Laird
|
||||
sum_w2_fixed = np.sum(w_fixed**2)
|
||||
numerator = chi2_fixed - dof
|
||||
denominator = sum_w_fixed - (sum_w2_fixed / sum_w_fixed)
|
||||
|
||||
# Tau^2 non può essere negativo. Se il Chi2 è minore del dof, Tau^2 = 0
|
||||
tau2 = max(0.0, numerator / denominator)
|
||||
|
||||
# --- STEP 3: Ricalcolo con Effetti Casuali (Random Effects) ---
|
||||
w_random = 1 / (ux**2 + tau2)
|
||||
sum_w_random = np.sum(w_random)
|
||||
|
||||
xm_random = np.sum(w_random * x) / sum_w_random
|
||||
uxm_random = 1 / np.sqrt(sum_w_random)
|
||||
|
||||
# Nota: Il p-value finale del Chi2 si riferisce alla consistenza iniziale dei dati
|
||||
p_value = stats.chi2.sf(chi2_fixed, dof)
|
||||
|
||||
return {
|
||||
"media": xm_random,
|
||||
"errore_media": uxm_random,
|
||||
"tau2": tau2,
|
||||
"tau": np.sqrt(tau2),
|
||||
"chi2_originale": chi2_fixed,
|
||||
"dof": dof,
|
||||
"p_value": p_value
|
||||
}
|
||||
|
||||
return
|
||||
|
||||
|
||||
@app.cell(hide_code=True)
|
||||
def _(mo):
|
||||
mo.md(r"""
|
||||
@@ -589,7 +679,7 @@ def _(mo):
|
||||
|
||||
@app.cell
|
||||
def _(np, plt, sns, stats):
|
||||
def plot_gauss(gaussiane):
|
||||
def plot_gauss(gaussiane, percorso=None, titolo="Confronto dei metodi"):
|
||||
# Creazione figura
|
||||
plt.figure(figsize=(12, 7))
|
||||
|
||||
@@ -600,10 +690,10 @@ def _(np, plt, sns, stats):
|
||||
for mu, sigma, _, _ in gaussiane:
|
||||
minimoLocale = mu - 4 * sigma
|
||||
massimoLocale = mu + 4 * sigma
|
||||
|
||||
|
||||
if minimoLocale < xMin:
|
||||
xMin = minimoLocale
|
||||
|
||||
|
||||
if massimoLocale > xMax:
|
||||
xMax = massimoLocale
|
||||
|
||||
@@ -613,7 +703,7 @@ def _(np, plt, sns, stats):
|
||||
for mu, sigma, colore, etichetta in gaussiane:
|
||||
y = stats.norm.pdf(x, mu, sigma)
|
||||
plt.plot(x, y, color=sns.color_palette()[colore], linewidth=1, label=etichetta)
|
||||
|
||||
|
||||
puntiLinee = [mu - sigma, mu, mu + sigma]
|
||||
for px in puntiLinee:
|
||||
py = stats.norm.pdf(px, mu, sigma)
|
||||
@@ -621,13 +711,17 @@ def _(np, plt, sns, stats):
|
||||
ymin=0, ymax=py,
|
||||
colors=sns.color_palette()[colore],
|
||||
linestyles='dashed', linewidth=1)
|
||||
|
||||
|
||||
# Dettagli estetici finali
|
||||
plt.ylim(bottom=0)
|
||||
plt.title('Confronto dati')
|
||||
plt.xlabel('k')
|
||||
plt.title(titolo)
|
||||
plt.xlabel('Focale Lente [$cm$]')
|
||||
plt.ylabel('Densità di probabilità [$cm^{-1}$]')
|
||||
plt.legend()
|
||||
|
||||
if percorso is not None:
|
||||
plt.savefig(str(percorso), dpi=200, bbox_inches="tight")
|
||||
|
||||
# Mostriamo il grafico
|
||||
plt.show()
|
||||
|
||||
|
||||
Reference in New Issue
Block a user