- sigma stimata solo con >= 8 campioni (altrimenti null); risultati.json senza NaN - adeguatezza = scala robusta dei residui (scala_residui_vs_rumore_max) - scarta_dopo_salto_n = 0, n_min_punti_parametri_fissi per i raffreddamenti a parametri fissi - salto del riscaldamento stimato dai modelli adiacenti, con ripiego sulla misura corretta Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
136 lines
7.5 KiB
Python
136 lines
7.5 KiB
Python
"""Analisi: fit dei modelli per ogni dominio e serie, scelta (docs/algoritmi.md §6) e stima dei salti del riscaldamento (§7)."""
|
||
|
||
import math
|
||
|
||
from .domini import sigma_rumore
|
||
from .fit import adatta
|
||
from .modelli import MODELLI, valuta
|
||
|
||
|
||
def valido(fit):
|
||
return (fit.get("convergenza") and not fit["parametri_a_limite"] and not fit["non_identificabile"]
|
||
and math.isfinite(fit["aicc"]))
|
||
|
||
|
||
def fit_scelto(vs):
|
||
"""Il fit scelto di una voce di serie, o None (indeterminato, riscaldamento)."""
|
||
return next((x for x in vs["fit"] if x["modello"] == vs["scelta"]), None)
|
||
|
||
|
||
def scegli(fit_list, sigma, cfg):
|
||
"""(fit scelto, motivo, adeguato); (None, motivo, False) se nessun fit è valido."""
|
||
f = cfg["fit"]
|
||
validi = [x for x in fit_list if valido(x)]
|
||
if not validi:
|
||
return None, "nessun modello convergente, identificabile e lontano dai limiti", False
|
||
noto = sigma is not None and math.isfinite(sigma) # con meno di 4 campioni σ non è stimabile
|
||
# adeguato = scala ROBUSTA dei residui vicina al rumore: come la perdita del fit, ignora gli scarti casuali (la RMSE no)
|
||
adeguati = [x for x in validi if noto and x["scala_residui"] <= f["scala_residui_vs_rumore_max"] * sigma]
|
||
pool = adeguati or validi
|
||
miglior = min(x["aicc"] for x in pool)
|
||
equivalenti = [x for x in pool if x["aicc"] <= miglior + f["soglia_delta_aicc"]]
|
||
scelto = min(equivalenti, key=lambda x: (x["n_parametri"], x["aicc"]))
|
||
motivo = (f"{len(validi)}/{len(fit_list)} modelli validi; "
|
||
+ ("modello adeguato" if adeguati else "nessun modello adeguato" if noto else "adeguatezza non valutabile (σ non stimabile)")
|
||
+ ("; tra i modelli equivalenti per AICc, quello con meno parametri" if len(equivalenti) > 1 else ""))
|
||
return scelto, motivo, bool(adeguati) if noto else None
|
||
|
||
|
||
def _gradino(t, y, d):
|
||
"""Descrizione del riscaldamento: salto misurato dal campione prima della salita all'ultimo campione di salita.
|
||
`_affina_salti` lo sostituisce con la stima dai modelli adiacenti quando possibile."""
|
||
i0, i1 = d["i0"], d["i1"]
|
||
prima, dopo = float(y[i0 - 1]), float(y[i1 - 1])
|
||
g = {"n": int(i1 - i0), "da": prima, "a": dopo, "delta": dopo - prima, "delta_misurato": dopo - prima,
|
||
"stima": "misura", "durata_s": float(t[i1 - 1] - t[i0 - 1])}
|
||
return {"sigma": None, "fit": [], "scelta": None, "gradino": g,
|
||
"motivo": f"salita rapida di {g['n']} campioni ({g['durata_s'] * 1000:.0f} ms), non fittata"}
|
||
|
||
|
||
def _serie(ty, yy, sigma, tipo, fissi, cfg):
|
||
voce = {"sigma": float(sigma) if math.isfinite(sigma) else None, "fit": [], "scelta": None}
|
||
chiave = "n_min_punti_parametri_fissi" if fissi else "n_min_punti" # con c e τ fissi basta stimare A: bastano meno campioni
|
||
minimo = cfg["domini"][chiave]
|
||
if len(yy) < minimo:
|
||
voce["motivo"] = f"{len(yy)} campioni usati nel fit, meno del minimo {minimo} (`{chiave}`)"
|
||
return voce
|
||
for m in cfg["fit"]["modelli_per_tipo"][tipo]:
|
||
nomi = MODELLI[m][1]
|
||
if len(nomi) - len([p for p in (fissi or {}) if p in nomi]) <= len(yy) - 2: # parametri liberi ≤ n − 2
|
||
voce["fit"].append(adatta(ty, yy, m, sigma, cfg, fissi))
|
||
scelto, voce["motivo"], voce["adeguato"] = scegli(voce["fit"], voce["sigma"], cfg)
|
||
voce["scelta"] = scelto["modello"] if scelto else None
|
||
return voce
|
||
|
||
|
||
def _valore(dominio, nome, t):
|
||
"""Valore del modello scelto di un dominio (per la serie `nome`) all'istante t; None se non ha un modello."""
|
||
fit = fit_scelto(dominio["serie"][nome])
|
||
if fit is None:
|
||
return None
|
||
return float(valuta(fit["modello"], t - dominio["t0_s"], list(fit["parametri"].values())))
|
||
|
||
|
||
def _finale(risultati):
|
||
"""Il risultato del raffreddamento finale (ultimo dominio di raffreddamento), o None."""
|
||
return next((r for r in reversed(risultati) if r["tipo"] == "raffreddamento"), None)
|
||
|
||
|
||
def _affina_salti(t, serie, risultati, domini):
|
||
"""Stima il salto di ogni riscaldamento all'istante del gradino (a metà tra l'ultimo campione prima della salita e il
|
||
primo in salita), perché la misura sui campioni lo sottostima (l'ultimo campione di salita è già in parte decaduto):
|
||
1. "modelli": differenza tra i modelli scelti per i due domini adiacenti, valutati in quell'istante;
|
||
2. "misura corretta": se manca un modello adiacente (tratto indeterminato) ma c'è il raffreddamento finale, la misura
|
||
con i due campioni riportati all'istante del gradino usando $c$ e τ del finale (stessa legge di decadimento);
|
||
3. "misura": altrimenti la misura sui campioni.
|
||
Una stima che si discosta più del 50 % dalla misura è scartata (guardia contro fit anomali)."""
|
||
finale = _finale(risultati)
|
||
for i, (r, d) in enumerate(zip(risultati, domini, strict=True)):
|
||
if r["tipo"] != "riscaldamento" or i == 0 or i + 1 >= len(risultati):
|
||
continue
|
||
ts = float(0.5 * (t[d["i0"] - 1] + t[d["i0"]]))
|
||
for nome, vs in r["serie"].items():
|
||
g = vs["gradino"]
|
||
prima, dopo = _valore(risultati[i - 1], nome, ts), _valore(risultati[i + 1], nome, ts)
|
||
stima = "modelli"
|
||
if prima is None or dopo is None:
|
||
fit = fit_scelto(finale["serie"][nome]) if finale else None
|
||
if fit is None or "tau" not in fit["parametri"]:
|
||
continue
|
||
c, tau, y = fit["parametri"]["c"], fit["parametri"]["tau"], serie[nome]
|
||
prima = c + (y[d["i0"] - 1] - c) * math.exp(-(ts - t[d["i0"] - 1]) / tau)
|
||
dopo = c + (y[d["i1"] - 1] - c) * math.exp((t[d["i1"] - 1] - ts) / tau)
|
||
stima = "misura corretta"
|
||
delta = float(dopo - prima)
|
||
if math.isfinite(delta) and abs(delta - g["delta_misurato"]) <= 0.5 * abs(g["delta_misurato"]):
|
||
g.update(da=float(prima), a=float(dopo), delta=delta, stima=stima, istante_s=ts)
|
||
|
||
|
||
def analizza_domini(t, serie, domini, cfg):
|
||
"""Un risultato per dominio. Il raffreddamento finale (ultimo) si fitta per primo, con parametri liberi; i
|
||
raffreddamenti intermedi (tra due gradini, corti) fissano i parametri di `raffreddamento_intermedio_fissa_da_finale`
|
||
ai suoi valori (stesse costanti fisiche). Il riscaldamento non si fitta: è descritto come salto."""
|
||
scarta = cfg["domini"]["scarta_dopo_salto_n"]
|
||
raffreddamenti = [i for i, d in enumerate(domini) if d["tipo"] == "raffreddamento"]
|
||
finale = raffreddamenti[-1] if raffreddamenti else None
|
||
risultati = [None] * len(domini)
|
||
fissi = {nome: {} for nome in serie}
|
||
for idx in sorted(range(len(domini)), key=lambda i: i != finale):
|
||
d = domini[idx]
|
||
i0, i1, t0 = d["i0"] + (scarta if d["tipo"] == "raffreddamento" else 0), d["i1"], t[d["i0"]]
|
||
voce = {"k": idx + 1, "tipo": d["tipo"], "t0_s": float(t0), "t1_s": float(t[i1 - 1]),
|
||
"n": int(i1 - d["i0"]), "n_fit": int(i1 - i0), "serie": {}}
|
||
for nome, y in serie.items():
|
||
if d["tipo"] == "riscaldamento":
|
||
voce["serie"][nome] = _gradino(t, y, d)
|
||
continue
|
||
yy = y[i0:i1]
|
||
usa = fissi[nome] if d["tipo"] == "raffreddamento" and idx != finale else None
|
||
vs = voce["serie"][nome] = _serie(t[i0:i1] - t0, yy, sigma_rumore(yy), d["tipo"], usa, cfg)
|
||
if idx == finale and vs["scelta"]:
|
||
par = next(x for x in vs["fit"] if x["modello"] == vs["scelta"])["parametri"]
|
||
fissi[nome] = {k: par[k] for k in cfg["fit"]["raffreddamento_intermedio_fissa_da_finale"] if k in par}
|
||
risultati[idx] = voce
|
||
_affina_salti(t, serie, risultati, domini)
|
||
return risultati
|