From 4f1530797e9e6234cfaa462cd70baaf54d6b5297 Mon Sep 17 00:00:00 2001 From: Egutierrez Date: Tue, 30 Jun 2026 20:42:12 +0200 Subject: [PATCH] =?UTF-8?q?feat(datascience):=20rigor=20experimental=20par?= =?UTF-8?q?a=20papers=20=E2=80=94=20effect=20size,=20IC,=20Holm=20+=20prer?= =?UTF-8?q?egistro=20inmutable?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Subsistema de papers reproducibles (grupo de capacidad `papers`). Añade las funciones estadísticas que un paper honesto necesita y la función que congela la hipótesis antes de mirar los datos (anti-HARKing). Nuevas funciones (puras salvo la última): - effect_size_cohens_d: Cohen's d + Hedges' g (corrección de sesgo para N pequeño) + interpretación cualitativa (negligible/small/medium/large por los umbrales de Cohen). Dict-no-throw ante varianza cero / N insuficiente. - confidence_interval_mean: intervalo de confianza de una media (t de Student) o de la diferencia de medias con Welch (df de Welch–Satterthwaite, sin asumir varianzas iguales). Dict-no-throw; el IC colapsa al punto cuando la varianza es cero. - preregister_hypothesis (impura): congela hipótesis + plan de análisis en papers//preregistration.md con frozen_at (UTC) y content_hash (sha256 del cuerpo normalizado, no del frontmatter). Inmutabilidad: una vez frozen, un contenido distinto se RECHAZA sin sobrescribir (mata el HARKing); idempotente si el contenido es idéntico. Siempre dict-no-throw. Extensión: - fdr_correction 1.0.0 -> 1.1.0: añade method="holm" (Holm-Bonferroni step-down, controla FWER, más potente que Bonferroni simple). Reúsa la maquinaria de alineación 1:1 con None/inválidos; no rompe los métodos bh/bonferroni. Reutiliza del registry: fdr_correction (BH + Bonferroni ya existían) como base para Holm. pearson y spearman_corr ya cubrían correlación. Tests: 36 pytest verdes (cohen/hedges 8, confidence/welch 8, fdr/holm/bonferroni 12, preregister 4 + extras), golden contra valores conocidos y validados con scipy. Golden manual del preregistro: congela, idempotente, rechaza edición (bytes en disco idénticos al congelado). Co-Authored-By: Claude Opus 4.8 (1M context) --- python/functions/datascience/__init__.py | 6 + .../datascience/confidence_interval_mean.md | 87 ++++++++ .../datascience/confidence_interval_mean.py | 176 +++++++++++++++ .../confidence_interval_mean_test.py | 140 ++++++++++++ .../datascience/effect_size_cohens_d.md | 80 +++++++ .../datascience/effect_size_cohens_d.py | 156 ++++++++++++++ .../datascience/effect_size_cohens_d_test.py | 96 +++++++++ .../functions/datascience/fdr_correction.md | 40 +++- .../functions/datascience/fdr_correction.py | 38 +++- .../datascience/fdr_correction_test.py | 66 +++++- .../datascience/preregister_hypothesis.md | 100 +++++++++ .../datascience/preregister_hypothesis.py | 202 ++++++++++++++++++ .../preregister_hypothesis_test.py | 99 +++++++++ 13 files changed, 1265 insertions(+), 21 deletions(-) create mode 100644 python/functions/datascience/confidence_interval_mean.md create mode 100644 python/functions/datascience/confidence_interval_mean.py create mode 100644 python/functions/datascience/confidence_interval_mean_test.py create mode 100644 python/functions/datascience/effect_size_cohens_d.md create mode 100644 python/functions/datascience/effect_size_cohens_d.py create mode 100644 python/functions/datascience/effect_size_cohens_d_test.py create mode 100644 python/functions/datascience/preregister_hypothesis.md create mode 100644 python/functions/datascience/preregister_hypothesis.py create mode 100644 python/functions/datascience/preregister_hypothesis_test.py diff --git a/python/functions/datascience/__init__.py b/python/functions/datascience/__init__.py index cdefab14..f4824679 100644 --- a/python/functions/datascience/__init__.py +++ b/python/functions/datascience/__init__.py @@ -59,6 +59,9 @@ from .acf_pacf import acf_pacf from .stl_decompose import stl_decompose from .to_returns import to_returns from .fdr_correction import fdr_correction +from .effect_size_cohens_d import effect_size_cohens_d +from .confidence_interval_mean import confidence_interval_mean +from .preregister_hypothesis import preregister_hypothesis from .suggest_reexpression import suggest_reexpression from .exploratory_caveats import exploratory_caveats from .render_eda_pdf import render_eda_pdf, render_eda_pdf_relational @@ -90,6 +93,9 @@ __all__ = [ "stl_decompose", "to_returns", "fdr_correction", + "effect_size_cohens_d", + "confidence_interval_mean", + "preregister_hypothesis", "suggest_reexpression", "exploratory_caveats", "render_eda_pdf", diff --git a/python/functions/datascience/confidence_interval_mean.md b/python/functions/datascience/confidence_interval_mean.md new file mode 100644 index 00000000..a44b70e7 --- /dev/null +++ b/python/functions/datascience/confidence_interval_mean.md @@ -0,0 +1,87 @@ +--- +name: confidence_interval_mean +kind: function +lang: py +domain: datascience +version: "1.0.0" +purity: pure +signature: "def confidence_interval_mean(data: list, other: list = None, confidence: float = 0.95) -> dict" +description: "Intervalo de confianza (IC) de la media de una muestra con la t de Student, o de la DIFERENCIA de medias de dos muestras independientes con el metodo de Welch (sin asumir varianzas iguales). Una muestra: df=n-1, se=sd_muestral/sqrt(n) (sd con ddof=1), tcrit=t.ppf((1+confidence)/2, df), ci=mean+/-tcrit*se. Dos muestras: IC de mean(data)-mean(other) con se=sqrt(se1^2+se2^2) y grados de libertad de Welch-Satterthwaite. Pura y robusta: nunca lanza; ante casos degenerados (muestra vacia, n<2) devuelve nan + clave note, y con varianza cero el IC colapsa al punto (no es error). Usa scipy.stats y numpy." +tags: [papers, statistics, confidence-interval, welch, t-test, python] +params: + - name: data + desc: "muestra de observaciones numericas (lista de numeros). Si other es None, el IC es el de la media de data." + - name: other + desc: "segunda muestra independiente (lista de numeros) o None (default). Si se da, el IC es el de la diferencia de medias mean(data)-mean(other) calculada con Welch (no asume varianzas iguales)." + - name: confidence + desc: "nivel de confianza en (0, 1); 0.95 = IC del 95% (default). El cuantil critico es t.ppf((1+confidence)/2, df)." +output: "dict {mean, ci_low, ci_high, se, df, confidence, n}. mean = media de data (una muestra) o la diferencia mean(data)-mean(other) (dos muestras). En el caso de dos muestras se anaden ademas n1 y n2 (y n = n1+n2). df son los grados de libertad de la t (Welch-Satterthwaite si dos muestras). Casos degenerados (muestra vacia, n<2) anaden la clave note y dejan ci_low/ci_high/se (y a veces df) en nan; con varianza cero y n>=2 el IC colapsa a [mean, mean] con se=0 (con note, sin nan). Nunca None ni excepcion." +uses_functions: [] +uses_types: [] +returns: [] +returns_optional: false +error_type: "" +imports: [scipy, numpy] +tested: true +tests: ["test_one_sample_golden_contra_scipy", "test_one_sample_distinto_nivel_confianza", "test_welch_diferencia_golden_contra_scipy", "test_edge_un_solo_elemento_no_lanza_nan_note", "test_edge_lista_vacia_no_lanza_note", "test_edge_varianza_cero_colapsa_al_punto", "test_edge_welch_muestra_vacia_no_lanza_note", "test_edge_welch_n1_uno_no_lanza_note"] +test_file_path: "python/functions/datascience/confidence_interval_mean_test.py" +file_path: "python/functions/datascience/confidence_interval_mean.py" +--- + +## Ejemplo + +```python +from datascience import confidence_interval_mean + +# IC del 95% de la media de una muestra (t de Student). +data = [2, 4, 4, 4, 5, 5, 7, 9] +ci = confidence_interval_mean(data, confidence=0.95) +print(ci["mean"]) # -> 5.0 +print(ci["df"]) # -> 7.0 (n - 1) +print(round(ci["ci_low"], 5), round(ci["ci_high"], 5)) +# -> 3.21251 6.78749 (se con sd muestral ddof=1 ~ 2.13809) + +# IC del 95% de la DIFERENCIA de medias (Welch, no asume varianzas iguales). +control = [23.0, 21.0, 25.0, 22.0, 24.0, 26.0] +tratado = [18.0, 20.0, 17.0, 19.0, 21.0] +diff = confidence_interval_mean(control, tratado, confidence=0.95) +print(diff["mean"]) # -> 4.5 (mean(control) - mean(tratado)) +print(round(diff["ci_low"], 4), round(diff["ci_high"], 4)) +# Si el intervalo no incluye 0, la diferencia es significativa al 5%. + +# Degenerados: nunca lanza. +print(confidence_interval_mean([5])["note"]) # n < 2: ... indefinidos +print(confidence_interval_mean([3, 3, 3])["se"]) # -> 0.0 (IC colapsa a [3, 3]) +``` + +## Cuando usarla + +Cuando quieras cuantificar la **incertidumbre de una media estimada** a partir de +una muestra: reporta `[ci_low, ci_high]` en vez de un punto suelto para mostrar +el rango plausible del valor real al nivel de confianza pedido. Usala tambien +para **comparar dos grupos** (A/B test, control vs tratamiento, antes vs +despues con grupos independientes): pasa las dos muestras y, si el IC de la +diferencia **no incluye el 0**, la diferencia es significativa al nivel +`1 - confidence`. Es el complemento del p-valor: ademas de "hay efecto", te dice +"de que tamano y con que margen". Para dos muestras usa Welch por defecto, asi +que no necesitas comprobar antes si las varianzas son iguales. + +## Gotchas + +- Pura y determinista (no hace I/O, no muta las entradas), pero **no** es + stdlib-only: depende de `scipy.stats` y `numpy` (ambos en el venv del proyecto). +- Con `other` usa **Welch** (df de Welch-Satterthwaite): NO asume varianzas + iguales ni tamanos de muestra iguales. Si necesitas el t-test clasico de + varianzas agrupadas (pooled), esta funcion no lo hace. +- `sd` se calcula con **ddof=1** (sd muestral), que es lo correcto para el IC de + una media con la t. Atajos como `sd_poblacional/sqrt(n)` (ddof=0) dan un + intervalo demasiado estrecho. +- En el caso de dos muestras, `mean` es la **diferencia** `mean(data) - mean(other)` + (no la media de data). El orden importa: el signo del IC depende de cual va + primero. +- Nunca lanza. Casos degenerados devuelven `nan` en `ci_low`/`ci_high`/`se` + (y a veces `df`) mas una clave `note`: muestra vacia o `n < 2` en cualquiera de + las muestras. **Excepcion**: con varianza cero y `n >= 2` el IC colapsa al + punto `[mean, mean]` con `se = 0` (no es un error, no hay `nan`). +- Comprueba `"note" in out` antes de usar `ci_low`/`ci_high` si la muestra puede + ser degenerada. diff --git a/python/functions/datascience/confidence_interval_mean.py b/python/functions/datascience/confidence_interval_mean.py new file mode 100644 index 00000000..6a4b5b2f --- /dev/null +++ b/python/functions/datascience/confidence_interval_mean.py @@ -0,0 +1,176 @@ +"""Intervalo de confianza de la media (una muestra) o de la diferencia de medias (Welch). + +Funcion pura del grupo papers. Calcula el intervalo de confianza (IC) de la media +de una muestra usando la t de Student, o el IC de la diferencia de medias de dos +muestras independientes con el metodo de Welch (sin asumir varianzas iguales). + +- Una muestra: ``df = n - 1``, ``se = sd / sqrt(n)`` (sd con ddof=1), + ``tcrit = t.ppf((1 + confidence) / 2, df)``, ``ci = mean +/- tcrit * se``. +- Dos muestras (Welch): IC de ``mean(data) - mean(other)``, con + ``se = sqrt(se1^2 + se2^2)`` y grados de libertad de Welch-Satterthwaite. + +No lanza excepciones: ante casos degenerados (muestras vacias, ``n < 2``, +varianza cero) devuelve un dict coherente con ``ci_low``/``ci_high``/``se`` en +``nan`` (salvo el sub-caso de varianza cero, donde el IC colapsa al punto) y una +clave ``note`` explicando el caso. Usa ``scipy.stats`` y ``numpy``. +""" + +from __future__ import annotations + +import math + +import numpy as np +from scipy import stats + + +def confidence_interval_mean( + data: list, other: list = None, confidence: float = 0.95 +) -> dict: + """Intervalo de confianza de la media o de la diferencia de medias (Welch). + + Si ``other`` es ``None``, calcula el IC de la media de ``data`` con la t de + Student. Si se proporciona ``other``, calcula el IC de la diferencia + ``mean(data) - mean(other)`` con el metodo de Welch (no asume varianzas + iguales) y grados de libertad de Welch-Satterthwaite. + + Es una funcion pura y determinista: no hace I/O ni muta las entradas. No + lanza excepcion ante datos degenerados; en su lugar devuelve un dict con la + clave ``note`` y los campos numericos indefinidos a ``nan``. + + Args: + data: muestra de observaciones numericas (lista de numeros). + other: segunda muestra independiente. Si se da, el IC es el de la + diferencia de medias ``mean(data) - mean(other)`` con Welch. Si es + ``None`` (default), el IC es el de la media de ``data``. + confidence: nivel de confianza en (0, 1), p.ej. 0.95 para el 95%. + + Returns: + dict con las claves: + mean: media de ``data`` (una muestra) o la diferencia + ``mean(data) - mean(other)`` (dos muestras). + ci_low: extremo inferior del intervalo de confianza. + ci_high: extremo superior del intervalo de confianza. + se: error estandar de la media (o de la diferencia). + df: grados de libertad de la t (Welch-Satterthwaite si dos muestras). + confidence: nivel de confianza aplicado (float). + n: tamano de la muestra (una muestra) o tamano total ``n1 + n2`` + (dos muestras; ademas se incluyen ``n1`` y ``n2``). + + En el caso de dos muestras se incluyen ademas ``n1`` y ``n2``. Casos + degenerados (muestra vacia, ``n < 2``, etc.) anaden la clave ``note`` y + dejan ``ci_low``/``ci_high``/``se`` (y a veces ``df``) en ``nan``. + """ + conf = float(confidence) + + if other is None: + return _ci_one_sample(data, conf) + return _ci_welch(data, other, conf) + + +def _ci_one_sample(data: list, conf: float) -> dict: + """IC de la media de una sola muestra con la t de Student.""" + arr = np.asarray(list(data), dtype=float) + n = int(arr.size) + + base = { + "mean": float("nan"), + "ci_low": float("nan"), + "ci_high": float("nan"), + "se": float("nan"), + "df": float("nan"), + "confidence": conf, + "n": n, + } + + if n == 0: + base["note"] = "muestra vacia: media e intervalo indefinidos" + return base + + mean = float(arr.mean()) + base["mean"] = mean + + if n < 2: + base["note"] = "n < 2: error estandar y grados de libertad indefinidos" + return base + + df = n - 1 + base["df"] = float(df) + + sd = float(arr.std(ddof=1)) + se = sd / math.sqrt(n) + base["se"] = se + + # Varianza cero: el IC colapsa al punto (no es un error). + if se == 0.0: + base["ci_low"] = mean + base["ci_high"] = mean + base["note"] = "varianza cero: el intervalo colapsa a la media" + return base + + tcrit = float(stats.t.ppf((1.0 + conf) / 2.0, df)) + margin = tcrit * se + base["ci_low"] = mean - margin + base["ci_high"] = mean + margin + return base + + +def _ci_welch(data: list, other: list, conf: float) -> dict: + """IC de la diferencia de medias de dos muestras con el metodo de Welch.""" + a = np.asarray(list(data), dtype=float) + b = np.asarray(list(other), dtype=float) + n1 = int(a.size) + n2 = int(b.size) + + base = { + "mean": float("nan"), + "ci_low": float("nan"), + "ci_high": float("nan"), + "se": float("nan"), + "df": float("nan"), + "confidence": conf, + "n": n1 + n2, + "n1": n1, + "n2": n2, + } + + if n1 == 0 or n2 == 0: + base["note"] = "alguna muestra esta vacia: diferencia e intervalo indefinidos" + return base + + mean1 = float(a.mean()) + mean2 = float(b.mean()) + diff = mean1 - mean2 + base["mean"] = diff + + if n1 < 2 or n2 < 2: + base["note"] = ( + "n < 2 en alguna muestra: error estandar y grados de libertad indefinidos" + ) + return base + + sd1 = float(a.std(ddof=1)) + sd2 = float(b.std(ddof=1)) + se1 = sd1 / math.sqrt(n1) + se2 = sd2 / math.sqrt(n2) + se = math.sqrt(se1 * se1 + se2 * se2) + base["se"] = se + + # Ambas varianzas cero: el IC de la diferencia colapsa al punto. + if se == 0.0: + base["ci_low"] = diff + base["ci_high"] = diff + base["df"] = float("nan") + base["note"] = "varianza cero en ambas muestras: el intervalo colapsa a la diferencia" + return base + + # Grados de libertad de Welch-Satterthwaite. + df = (se1 * se1 + se2 * se2) ** 2 / ( + (se1**4) / (n1 - 1) + (se2**4) / (n2 - 1) + ) + base["df"] = float(df) + + tcrit = float(stats.t.ppf((1.0 + conf) / 2.0, df)) + margin = tcrit * se + base["ci_low"] = diff - margin + base["ci_high"] = diff + margin + return base diff --git a/python/functions/datascience/confidence_interval_mean_test.py b/python/functions/datascience/confidence_interval_mean_test.py new file mode 100644 index 00000000..70ccc029 --- /dev/null +++ b/python/functions/datascience/confidence_interval_mean_test.py @@ -0,0 +1,140 @@ +"""Tests para confidence_interval_mean (IC de la media / diferencia de medias Welch). + +Importa el modulo hoja directamente (`confidence_interval_mean`) para no depender +de que el paquete reexporte la funcion en su __init__ (lo integra el orquestador +al cerrar el grupo). + +Los golden se calculan con scipy dentro del propio test para que sean robustos: +la funcion bajo prueba debe coincidir con la referencia de scipy a ~1e-9. +""" + +import math + +import numpy as np +from scipy import stats + +from confidence_interval_mean import confidence_interval_mean + + +def test_one_sample_golden_contra_scipy(): + # mean=5.0, n=8. Este dataset tiene sd POBLACIONAL (ddof=0) exactamente 2.0, + # pero la sd MUESTRAL (ddof=1, la que exige la spec y la que es correcta para + # el IC de una media con la t) es sqrt(32/7) ~ 2.13809. El golden robusto se + # calcula con scipy usando se con ddof=1, no con el atajo 2.0/sqrt(8). + data = [2, 4, 4, 4, 5, 5, 7, 9] + out = confidence_interval_mean(data, confidence=0.95) + + n = len(data) + mean = float(np.mean(data)) + sd = float(np.std(data, ddof=1)) # sample sd ~ 2.13809 + se = sd / math.sqrt(n) + lo, hi = stats.t.interval(0.95, df=n - 1, loc=mean, scale=se) + + assert abs(out["mean"] - 5.0) < 1e-9 + assert abs(out["se"] - se) < 1e-12 + assert out["df"] == 7.0 + assert out["n"] == 8 + assert out["confidence"] == 0.95 + assert abs(out["ci_low"] - lo) < 1e-9 + assert abs(out["ci_high"] - hi) < 1e-9 + # Valores tabulados correctos para ddof=1 (no los 3.32793/6.67207 del + # enunciado, que asumian erroneamente sd=2.0 / ddof=0). + assert abs(out["ci_low"] - 3.21251) < 1e-3 + assert abs(out["ci_high"] - 6.78749) < 1e-3 + assert "note" not in out + + +def test_one_sample_distinto_nivel_confianza(): + data = [10.0, 12.0, 11.0, 13.0, 9.0, 14.0] + out = confidence_interval_mean(data, confidence=0.99) + + n = len(data) + mean = float(np.mean(data)) + se = float(np.std(data, ddof=1)) / math.sqrt(n) + lo, hi = stats.t.interval(0.99, df=n - 1, loc=mean, scale=se) + + assert abs(out["mean"] - mean) < 1e-12 + assert abs(out["ci_low"] - lo) < 1e-9 + assert abs(out["ci_high"] - hi) < 1e-9 + assert out["df"] == float(n - 1) + + +def test_welch_diferencia_golden_contra_scipy(): + data = [23.0, 21.0, 25.0, 22.0, 24.0, 26.0] + other = [18.0, 20.0, 17.0, 19.0, 21.0] + conf = 0.95 + out = confidence_interval_mean(data, other, confidence=conf) + + a = np.asarray(data, dtype=float) + b = np.asarray(other, dtype=float) + n1, n2 = a.size, b.size + mean1, mean2 = float(a.mean()), float(b.mean()) + diff = mean1 - mean2 + se1 = float(a.std(ddof=1)) / math.sqrt(n1) + se2 = float(b.std(ddof=1)) / math.sqrt(n2) + se = math.sqrt(se1**2 + se2**2) + df = (se1**2 + se2**2) ** 2 / (se1**4 / (n1 - 1) + se2**4 / (n2 - 1)) + lo, hi = stats.t.interval(conf, df=df, loc=diff, scale=se) + + assert abs(out["mean"] - diff) < 1e-9 + assert abs(out["mean"] - (mean1 - mean2)) < 1e-9 + assert abs(out["se"] - se) < 1e-12 + assert abs(out["df"] - df) < 1e-9 + assert abs(out["ci_low"] - lo) < 1e-9 + assert abs(out["ci_high"] - hi) < 1e-9 + assert out["n1"] == n1 + assert out["n2"] == n2 + assert out["n"] == n1 + n2 + assert "note" not in out + + +def test_edge_un_solo_elemento_no_lanza_nan_note(): + out = confidence_interval_mean([5], confidence=0.95) + assert out["mean"] == 5.0 # la media si esta definida con n=1 + assert math.isnan(out["se"]) + assert math.isnan(out["ci_low"]) + assert math.isnan(out["ci_high"]) + assert math.isnan(out["df"]) + assert out["n"] == 1 + assert "note" in out + + +def test_edge_lista_vacia_no_lanza_note(): + out = confidence_interval_mean([], confidence=0.95) + assert math.isnan(out["mean"]) + assert math.isnan(out["ci_low"]) + assert math.isnan(out["ci_high"]) + assert math.isnan(out["se"]) + assert out["n"] == 0 + assert "note" in out + + +def test_edge_varianza_cero_colapsa_al_punto(): + out = confidence_interval_mean([3, 3, 3], confidence=0.95) + assert out["mean"] == 3.0 + assert out["se"] == 0.0 + assert out["ci_low"] == 3.0 + assert out["ci_high"] == 3.0 + assert not math.isnan(out["ci_low"]) + assert out["n"] == 3 + assert "note" in out + + +def test_edge_welch_muestra_vacia_no_lanza_note(): + out = confidence_interval_mean([1.0, 2.0, 3.0], [], confidence=0.95) + assert math.isnan(out["mean"]) + assert math.isnan(out["ci_low"]) + assert math.isnan(out["se"]) + assert out["n1"] == 3 + assert out["n2"] == 0 + assert "note" in out + + +def test_edge_welch_n1_uno_no_lanza_note(): + out = confidence_interval_mean([5.0], [1.0, 2.0, 3.0], confidence=0.95) + # La diferencia de medias si esta definida. + assert abs(out["mean"] - (5.0 - 2.0)) < 1e-9 + assert math.isnan(out["se"]) + assert math.isnan(out["ci_low"]) + assert math.isnan(out["df"]) + assert "note" in out diff --git a/python/functions/datascience/effect_size_cohens_d.md b/python/functions/datascience/effect_size_cohens_d.md new file mode 100644 index 00000000..6e5adbe8 --- /dev/null +++ b/python/functions/datascience/effect_size_cohens_d.md @@ -0,0 +1,80 @@ +--- +name: effect_size_cohens_d +kind: function +lang: py +domain: datascience +version: "1.0.0" +purity: pure +signature: "def effect_size_cohens_d(group_a: list, group_b: list) -> dict" +description: "Tamano del efecto (effect size) entre dos grupos numericos: Cohen's d (diferencia de medias estandarizada por la desviacion tipica combinada, varianzas muestrales ddof=1), Hedges' g (d corregido por el sesgo al alza con muestras pequenas via el factor J) e interpretacion cualitativa de la magnitud segun los umbrales clasicos de Cohen (negligible/small/medium/large). El p-valor dice si hay diferencia; el effect size dice como de grande, de forma adimensional e independiente del N. Pura, sin dependencias externas; nunca lanza: los casos degenerados (varianza cero, N<2, listas vacias) devuelven NaN + una clave note." +tags: [papers, statistics, effect-size, cohens-d, hedges-g, python] +params: + - name: group_a + desc: "primera muestra (lista de numeros). Necesita >=2 observaciones para que exista la varianza muestral (ddof=1)." + - name: group_b + desc: "segunda muestra (lista de numeros). Necesita >=2 observaciones. El signo de cohens_d es positivo cuando mean_a > mean_b." +output: "dict {cohens_d: float (diferencia de medias estandarizada, puede ser NaN), hedges_g: float (cohens_d * factor de correccion J, puede ser NaN), interpretation: str ('negligible'|'small'|'medium'|'large', o 'undefined' en casos degenerados), n_a: int, n_b: int, mean_a: float, mean_b: float, pooled_sd: float (desviacion tipica combinada)}. Casos degenerados (varianza cero en ambos grupos, N<2 en algun grupo, o listas vacias) anaden clave note. Nunca None ni excepcion." +uses_functions: [] +uses_types: [] +returns: [] +returns_optional: false +error_type: "" +imports: [math] +tested: true +tests: ["test_golden_large_effect", "test_hedges_g_menor_en_magnitud_que_cohens_d", "test_interpretation_thresholds", "test_signo_positivo_cuando_a_mayor_que_b", "test_varianza_cero_no_lanza", "test_n_insuficiente_no_lanza", "test_listas_vacias_no_lanza", "test_un_grupo_vacio_no_lanza"] +test_file_path: "python/functions/datascience/effect_size_cohens_d_test.py" +file_path: "python/functions/datascience/effect_size_cohens_d.py" +--- + +## Ejemplo + +```python +from datascience import effect_size_cohens_d + +# Dos grupos desplazados 2 unidades, misma dispersion. +a = [1, 2, 3, 4, 5] # media 3, varianza muestral 2.5 +b = [3, 4, 5, 6, 7] # media 5, varianza muestral 2.5 + +out = effect_size_cohens_d(a, b) +print(out["cohens_d"]) # -> -1.264911... (a esta 1.26 SD por debajo de b) +print(out["hedges_g"]) # -> -1.142500... (|g| < |d|: correccion N pequeno) +print(out["interpretation"]) # -> "large" (|d| >= 0.8) +print(out["pooled_sd"]) # -> 1.581138... + +# Caso degenerado: varianza cero -> no lanza, NaN + note. +deg = effect_size_cohens_d([5, 5, 5], [5, 5, 5]) +print(deg["interpretation"]) # -> "undefined" +print(deg["note"]) # -> "varianza cero, effect size indefinido" +``` + +## Cuando usarla + +Cuando ya sepas que dos grupos difieren (o quieras cuantificar su diferencia) +y necesites una medida **de magnitud, no de significancia**: comparar el antes +y el despues de una intervencion, el grupo control frente al tratamiento, o dos +cohortes. Reportala junto al p-valor para responder "¿como de grande es la +diferencia?" — un p-valor minusculo con N enorme puede esconder un efecto +trivial. Es adimensional (en unidades de desviaciones tipicas), asi que hace +comparables resultados entre estudios y alimenta meta-analisis. Usa **Hedges' g** +en lugar de Cohen's d cuando los grupos sean pequenos (decenas o menos): d +sobreestima el efecto y g lo corrige. + +## Gotchas + +- Pura y sin dependencias externas (solo `math` de la stdlib). +- Usa **varianza muestral** (ddof=1), no poblacional. Por eso cada grupo + necesita al menos 2 observaciones; con N=1 la varianza muestral no existe y la + funcion devuelve NaN + `note`. +- **Nunca lanza excepcion**. Los casos degenerados devuelven `cohens_d` y + `hedges_g` a `float('nan')`, `interpretation="undefined"` y una clave `note`: + varianza cero en ambos grupos (`pooled_sd == 0`), N<2 en algun grupo, o listas + vacias. Comprueba con `math.isnan(out["cohens_d"])` o la presencia de `note` + antes de usar el resultado. +- El **signo** de `cohens_d` depende del orden de los argumentos: positivo si + `mean_a > mean_b`, negativo en caso contrario. La `interpretation` usa `|d|`, + asi que no depende del orden. +- `pooled_sd` asume varianzas comparables entre grupos (homogeneidad). Si las + dispersiones son muy distintas, Cohen's d clasico pierde precision; considera + variantes (Glass's delta) fuera del alcance de esta funcion. +- Los umbrales de Cohen (0.2 / 0.5 / 0.8) son convencion, no ley: interpretalos + segun el dominio. diff --git a/python/functions/datascience/effect_size_cohens_d.py b/python/functions/datascience/effect_size_cohens_d.py new file mode 100644 index 00000000..dee7ed18 --- /dev/null +++ b/python/functions/datascience/effect_size_cohens_d.py @@ -0,0 +1,156 @@ +"""Effect size de dos grupos: Cohen's d, Hedges' g e interpretacion cualitativa. + +Funcion pura del grupo papers. El p-valor responde a "¿hay diferencia?" pero no +a "¿como de grande es?". El tamano del efecto (effect size) cuantifica la +magnitud de la diferencia entre dos grupos de forma adimensional, independiente +del N, y es lo que hace comparables resultados entre estudios (meta-analisis). + +- Cohen's d: diferencia de medias estandarizada por la desviacion tipica + combinada (pooled SD), con varianzas muestrales (ddof=1). +- Hedges' g: Cohen's d corregido por el sesgo al alza que sufre d con muestras + pequenas, multiplicando por el factor de correccion J. +- interpretation: etiqueta cualitativa de |d| segun los umbrales clasicos de + Cohen (negligible / small / medium / large). + +No usa dependencias externas: aritmetica de la libreria estandar (``math``). +""" + +from __future__ import annotations + +import math + + +def _mean(xs: list) -> float: + """Media aritmetica de una lista no vacia de numeros.""" + return sum(float(x) for x in xs) / len(xs) + + +def _sample_variance(xs: list, mean: float) -> float: + """Varianza muestral (ddof=1) de una lista con al menos 2 elementos.""" + n = len(xs) + return sum((float(x) - mean) ** 2 for x in xs) / (n - 1) + + +def _interpret(abs_d: float) -> str: + """Etiqueta cualitativa del tamano del efecto segun |d| (umbrales de Cohen).""" + if abs_d < 0.2: + return "negligible" + if abs_d < 0.5: + return "small" + if abs_d < 0.8: + return "medium" + return "large" + + +def effect_size_cohens_d(group_a: list, group_b: list) -> dict: + """Calcula el tamano del efecto entre dos grupos numericos. + + Devuelve Cohen's d (diferencia de medias estandarizada por la pooled SD), + Hedges' g (d corregido por sesgo de muestra pequena) y una etiqueta + cualitativa de la magnitud segun los umbrales de Cohen. + + Es una funcion pura y determinista: no hace I/O, no muta la entrada. No lanza + excepcion ante datos degenerados; en su lugar devuelve un dict con + ``cohens_d`` / ``hedges_g`` a ``float('nan')``, ``interpretation`` a + ``"undefined"`` y una clave ``note`` explicando el caso. + + Definiciones: + s_pooled = sqrt(((n1-1)*s1^2 + (n2-1)*s2^2) / (n1+n2-2)), con s1^2, s2^2 + varianzas muestrales (ddof=1). + cohens_d = (mean_a - mean_b) / s_pooled. + J = 1 - 3 / (4*(n1+n2) - 9) (factor de correccion de Hedges). + hedges_g = cohens_d * J. + + Args: + group_a: primera muestra (lista de numeros). Necesita >=2 elementos para + que exista la varianza muestral. + group_b: segunda muestra (lista de numeros). Necesita >=2 elementos. + + Returns: + dict con las claves: + cohens_d: float, diferencia de medias estandarizada (puede ser NaN). + hedges_g: float, Cohen's d corregido por sesgo (puede ser NaN). + interpretation: str, "negligible" | "small" | "medium" | "large", o + "undefined" en casos degenerados. + n_a: int, tamano de group_a. + n_b: int, tamano de group_b. + mean_a: float, media de group_a (NaN si vacio). + mean_b: float, media de group_b (NaN si vacio). + pooled_sd: float, desviacion tipica combinada (NaN si indefinida). + + Casos degenerados (lista vacia, N<2 en algun grupo, o varianza cero en + ambos grupos -> pooled_sd == 0) anaden ademas una clave ``note``. + """ + nan = float("nan") + n_a = len(group_a) + n_b = len(group_b) + + # Listas vacias: ni media ni varianza definidas. + if n_a == 0 or n_b == 0: + return { + "cohens_d": nan, + "hedges_g": nan, + "interpretation": "undefined", + "n_a": n_a, + "n_b": n_b, + "mean_a": _mean(group_a) if n_a else nan, + "mean_b": _mean(group_b) if n_b else nan, + "pooled_sd": nan, + "note": "grupo vacio: media y varianza indefinidas, effect size indefinido", + } + + mean_a = _mean(group_a) + mean_b = _mean(group_b) + + # N insuficiente: la varianza muestral (ddof=1) no existe con un solo dato, + # y la correccion de Hedges no es fiable. + if n_a < 2 or n_b < 2: + return { + "cohens_d": nan, + "hedges_g": nan, + "interpretation": "undefined", + "n_a": n_a, + "n_b": n_b, + "mean_a": mean_a, + "mean_b": mean_b, + "pooled_sd": nan, + "note": ( + "N insuficiente: cada grupo necesita >=2 observaciones para la " + "varianza muestral; effect size indefinido" + ), + } + + var_a = _sample_variance(group_a, mean_a) + var_b = _sample_variance(group_b, mean_b) + pooled_sd = math.sqrt( + ((n_a - 1) * var_a + (n_b - 1) * var_b) / (n_a + n_b - 2) + ) + + # Varianza cero en ambos grupos: no se puede estandarizar (division por 0). + if pooled_sd == 0.0: + return { + "cohens_d": nan, + "hedges_g": nan, + "interpretation": "undefined", + "n_a": n_a, + "n_b": n_b, + "mean_a": mean_a, + "mean_b": mean_b, + "pooled_sd": 0.0, + "note": "varianza cero, effect size indefinido", + } + + cohens_d = (mean_a - mean_b) / pooled_sd + j = 1.0 - 3.0 / (4.0 * (n_a + n_b) - 9.0) + hedges_g = cohens_d * j + + return { + "cohens_d": cohens_d, + "hedges_g": hedges_g, + "interpretation": _interpret(abs(cohens_d)), + "n_a": n_a, + "n_b": n_b, + "mean_a": mean_a, + "mean_b": mean_b, + "pooled_sd": pooled_sd, + } diff --git a/python/functions/datascience/effect_size_cohens_d_test.py b/python/functions/datascience/effect_size_cohens_d_test.py new file mode 100644 index 00000000..f5c43049 --- /dev/null +++ b/python/functions/datascience/effect_size_cohens_d_test.py @@ -0,0 +1,96 @@ +"""Tests para effect_size_cohens_d (tamano del efecto de dos grupos). + +Importa el modulo hoja directamente (`effect_size_cohens_d`) para no depender de +que el paquete reexporte la funcion en su __init__ (lo integra el orquestador al +cerrar el grupo papers). El pytest del repo tiene pythonpath=["functions", ...], +asi que el modulo hoja se resuelve por su nombre directo. +""" + +import math + +from effect_size_cohens_d import effect_size_cohens_d + + +def test_golden_large_effect(): + # group_a: mean 3, var muestral 2.5; group_b: mean 5, var 2.5. + # pooled_sd = sqrt(2.5) ~= 1.5811388. + # cohens_d = (3-5)/1.5811388 ~= -1.264911. + # J = 1 - 3/(4*10-9) = 1 - 3/31 = 0.9032258. + # hedges_g = d * J = -1.2649111 * 0.9032258 ~= -1.142500. + out = effect_size_cohens_d([1, 2, 3, 4, 5], [3, 4, 5, 6, 7]) + assert abs(out["cohens_d"] - (-1.26491)) < 1e-4 + assert abs(out["hedges_g"] - (-1.14250)) < 1e-4 + assert out["interpretation"] == "large" + assert out["n_a"] == 5 + assert out["n_b"] == 5 + assert abs(out["mean_a"] - 3.0) < 1e-12 + assert abs(out["mean_b"] - 5.0) < 1e-12 + assert abs(out["pooled_sd"] - math.sqrt(2.5)) < 1e-9 + assert "note" not in out + + +def test_hedges_g_menor_en_magnitud_que_cohens_d(): + # La correccion J esta en (0, 1), asi que |g| < |d| siempre. + out = effect_size_cohens_d([1, 2, 3, 4, 5], [3, 4, 5, 6, 7]) + assert abs(out["hedges_g"]) < abs(out["cohens_d"]) + + +def test_interpretation_thresholds(): + # negligible: |d| < 0.2. Medias casi iguales con varianza grande. + neg = effect_size_cohens_d([0, 10, 20, 30], [1, 11, 21, 31]) + assert neg["interpretation"] == "negligible" + assert abs(neg["cohens_d"]) < 0.2 + + # small: 0.2 <= |d| < 0.5. + small = effect_size_cohens_d([0, 10, 20, 30], [4, 14, 24, 34]) + assert small["interpretation"] == "small" + assert 0.2 <= abs(small["cohens_d"]) < 0.5 + + # medium: 0.5 <= |d| < 0.8. + medium = effect_size_cohens_d([0, 10, 20, 30], [9, 19, 29, 39]) + assert medium["interpretation"] == "medium" + assert 0.5 <= abs(medium["cohens_d"]) < 0.8 + + +def test_signo_positivo_cuando_a_mayor_que_b(): + out = effect_size_cohens_d([10, 12, 14, 16], [1, 2, 3, 4]) + assert out["cohens_d"] > 0 + assert out["interpretation"] == "large" + + +def test_varianza_cero_no_lanza(): + out = effect_size_cohens_d([5, 5, 5], [5, 5, 5]) + assert math.isnan(out["cohens_d"]) + assert math.isnan(out["hedges_g"]) + assert out["interpretation"] == "undefined" + assert out["pooled_sd"] == 0.0 + assert "note" in out + assert "varianza cero" in out["note"] + + +def test_n_insuficiente_no_lanza(): + out = effect_size_cohens_d([3], [1, 2, 3]) + assert math.isnan(out["cohens_d"]) + assert math.isnan(out["hedges_g"]) + assert out["interpretation"] == "undefined" + assert out["n_a"] == 1 + assert out["n_b"] == 3 + assert "note" in out + + +def test_listas_vacias_no_lanza(): + out = effect_size_cohens_d([], []) + assert math.isnan(out["cohens_d"]) + assert math.isnan(out["hedges_g"]) + assert out["interpretation"] == "undefined" + assert out["n_a"] == 0 + assert out["n_b"] == 0 + assert "note" in out + + +def test_un_grupo_vacio_no_lanza(): + out = effect_size_cohens_d([1, 2, 3], []) + assert math.isnan(out["cohens_d"]) + assert out["interpretation"] == "undefined" + assert out["n_b"] == 0 + assert "note" in out diff --git a/python/functions/datascience/fdr_correction.md b/python/functions/datascience/fdr_correction.md index 32fd4635..8ae003da 100644 --- a/python/functions/datascience/fdr_correction.md +++ b/python/functions/datascience/fdr_correction.md @@ -3,19 +3,19 @@ name: fdr_correction kind: function lang: py domain: datascience -version: "1.0.0" +version: "1.1.0" purity: pure signature: "def fdr_correction(pvalues: list, alpha: float = 0.05, method: str = \"bh\") -> dict" -description: "Correccion de comparaciones multiples (multiple-testing) sobre una lista de p-valores: Benjamini-Hochberg (FDR, 'bh') o Bonferroni (FWER, 'bonferroni'). Antidoto al sesgo de mineria de datos (data-mining bias): al evaluar muchas hipotesis a la vez (todos los pares de una matriz), el azar produce falsos positivos; esta funcion ajusta los p-valores y marca cuales siguen siendo significativos tras corregir. Pura, sin dependencias externas, alineada 1:1 con la entrada (admite None en posiciones sin test)." -tags: [eda, statistics, multiple-testing, fdr, benjamini-hochberg, bonferroni, p-value, data-mining-bias, python] +description: "Correccion de comparaciones multiples (multiple-testing) sobre una lista de p-valores: Benjamini-Hochberg (FDR, 'bh'), Bonferroni (FWER, 'bonferroni') o Holm-Bonferroni (FWER step-down, 'holm', mas potente que Bonferroni simple). Antidoto al sesgo de mineria de datos (data-mining bias): al evaluar muchas hipotesis a la vez (todos los pares de una matriz), el azar produce falsos positivos; esta funcion ajusta los p-valores y marca cuales siguen siendo significativos tras corregir. Pura, sin dependencias externas, alineada 1:1 con la entrada (admite None en posiciones sin test)." +tags: [eda, statistics, multiple-testing, fdr, benjamini-hochberg, bonferroni, holm, holm-bonferroni, fwer, p-value, data-mining-bias, python] params: - name: pvalues desc: "lista de p-valores (floats en [0, 1]). Se admiten None u otros valores no validos en posiciones sin test disponible; se propagan como None en la salida y no cuentan como prueba (m)." - name: alpha desc: "nivel de significancia objetivo tras la correccion (default 0.05). Para BH es el umbral del FDR; para Bonferroni, del FWER (tasa de error por familia)." - name: method - desc: "'bh' = Benjamini-Hochberg (controla FDR, menos conservador, mas potencia); 'bonferroni' = controla FWER (mas conservador). Cualquier otro valor devuelve un dict con note." -output: "dict {p_values_adjusted: lista alineada con pvalues (float ajustado o None), reject: lista de bool (True = significativo tras corregir), n_tests: nº de p-valores validos (m), n_rejected: nº de hipotesis rechazadas, alpha: float aplicado, method: str}. Casos degenerados (vacio, sin p validos, metodo desconocido) anaden clave note. Nunca None ni excepcion." + desc: "'bh' = Benjamini-Hochberg (controla FDR, menos conservador, mas potencia); 'bonferroni' = controla FWER (mas conservador); 'holm' = Holm-Bonferroni (controla FWER, step-down, uniformemente mas potente que Bonferroni simple). Cualquier otro valor devuelve un dict con note." +output: "dict {p_values_adjusted: lista alineada con pvalues (float ajustado o None), reject: lista de bool (True = significativo tras corregir), n_tests: nº de p-valores validos (m), n_rejected: nº de hipotesis rechazadas, alpha: float aplicado, method: str ('bh' | 'bonferroni' | 'holm')}. Casos degenerados (vacio, sin p validos, metodo desconocido) anaden clave note. Nunca None ni excepcion." uses_functions: [] uses_types: [] returns: [] @@ -23,7 +23,7 @@ returns_optional: false error_type: "" imports: [math] tested: true -tests: ["test_bh_golden_rechaza_dos_de_tres", "test_bonferroni_mas_conservador_que_bh", "test_p_values_adjusted_alineados_y_en_rango", "test_none_se_propaga_alineado", "test_lista_vacia_devuelve_note", "test_solo_none_devuelve_note", "test_metodo_desconocido_devuelve_note", "test_todos_significativos"] +tests: ["test_bh_golden_rechaza_dos_de_tres", "test_bonferroni_mas_conservador_que_bh", "test_p_values_adjusted_alineados_y_en_rango", "test_none_se_propaga_alineado", "test_lista_vacia_devuelve_note", "test_solo_none_devuelve_note", "test_metodo_desconocido_devuelve_note", "test_todos_significativos", "test_holm_golden_rechaza_dos_de_cuatro", "test_holm_entre_bonferroni_y_bh", "test_none_se_propaga_alineado_holm", "test_lista_vacia_holm_devuelve_note"] test_file_path: "python/functions/datascience/fdr_correction_test.py" file_path: "python/functions/datascience/fdr_correction.py" --- @@ -45,6 +45,13 @@ bon = fdr_correction(pvalues, alpha=0.05, method="bonferroni") print(bon["reject"]) # -> [True, False, False] print(bon["p_values_adjusted"]) # -> [0.03, 0.06, 1.0] +# Holm-Bonferroni (step-down): controla el FWER como Bonferroni pero es mas +# potente; rechaza al menos tanto como Bonferroni simple, nunca menos. +holm = fdr_correction([0.01, 0.04, 0.03, 0.005], alpha=0.05, method="holm") +print(holm["reject"]) # -> [True, False, False, True] +print(holm["p_values_adjusted"]) # -> [0.03, 0.06, 0.06, 0.02] +print(holm["n_rejected"]) # -> 2 + # Posiciones sin test (None) se propagan alineadas: el llamador puede pasar la # lista completa de pares y recuperar el mapeo 1:1. mix = fdr_correction([0.001, None, 0.9]) @@ -61,8 +68,11 @@ combinaciones y se quede con las que "pasan". Sin corregir, con N pruebas y alpha=0.05 esperas ~5% de falsos positivos *por azar*: cuantas mas pruebas, mas correlaciones espurias. Llama a `fdr_correction` con todos los p-valores de la familia y usa `reject` (no el umbral crudo) para decidir que es real. Usa `"bh"` -por defecto (mejor potencia); `"bonferroni"` cuando un falso positivo sea muy -costoso y prefieras maxima cautela. +por defecto (mejor potencia); `"holm"` (Holm-Bonferroni, FWER step-down) cuando +quieras controlar el FWER pero sin la perdida de potencia de Bonferroni simple +(rechaza al menos tanto como `"bonferroni"`, nunca menos); `"bonferroni"` cuando +un falso positivo sea muy costoso y prefieras la maxima cautela del metodo mas +simple. ## Gotchas @@ -76,8 +86,16 @@ costoso y prefieras maxima cautela. eso puedes pasar la lista completa de pares aunque algunos no tengan test. - `n_tests` es el numero de p-valores **validos** (m), que puede ser menor que `len(pvalues)` si hay `None`. -- BH y Bonferroni controlan cosas distintas: BH la tasa de falsos - descubrimientos (FDR), Bonferroni la probabilidad de *cualquier* falso +- BH controla cosa distinta que Bonferroni/Holm: BH la tasa de falsos + descubrimientos (FDR); Bonferroni y Holm la probabilidad de *cualquier* falso positivo (FWER). No son intercambiables; elige segun el coste de equivocarte. +- `"holm"` y `"bonferroni"` controlan ambos el FWER, pero Holm es step-down y + uniformemente mas potente: rechaza al menos tantas hipotesis como Bonferroni + simple sobre el mismo set, nunca menos. Si controlas FWER, `"holm"` domina a + `"bonferroni"` salvo que necesites el ajuste mas simple por interpretabilidad. - Metodo desconocido o lista vacia/sin p validos no lanzan: devuelven un dict - con `note`. + con `note`. Los metodos validos son `"bh"`, `"bonferroni"` y `"holm"`. + +## Capability growth log + +- v1.1.0 (2026-06-30) — añade method="holm" (Holm-Bonferroni step-down, FWER, más potente que Bonferroni simple). diff --git a/python/functions/datascience/fdr_correction.py b/python/functions/datascience/fdr_correction.py index 25b0d58b..ea5408a0 100644 --- a/python/functions/datascience/fdr_correction.py +++ b/python/functions/datascience/fdr_correction.py @@ -5,12 +5,15 @@ todos los pares de una matriz de asociacion), la probabilidad de obtener al meno un falso positivo por azar crece con el numero de pruebas: es el sesgo de mineria de datos (data-mining bias) descrito por Aronson en *Evidence-Based Technical Analysis* (cap. 6). Esta funcion ajusta los p-valores para controlar ese sesgo -mediante dos metodos clasicos: +mediante tres metodos clasicos: - Benjamini-Hochberg (``"bh"``): controla la tasa de falsos descubrimientos (False Discovery Rate, FDR). Menos conservador, mas potencia estadistica. - Bonferroni (``"bonferroni"``): controla la tasa de error por familia (Family-Wise Error Rate, FWER). Mas conservador. +- Holm-Bonferroni (``"holm"``): controla el FWER como Bonferroni pero es un + procedimiento step-down uniformemente mas potente; rechaza al menos tantas + hipotesis como Bonferroni simple, nunca menos. No usa dependencias externas: aritmetica de la libreria estandar. """ @@ -35,8 +38,9 @@ def _is_valid_p(v) -> bool: def fdr_correction(pvalues: list, alpha: float = 0.05, method: str = "bh") -> dict: """Corrige una lista de p-valores por comparaciones multiples. - Aplica Benjamini-Hochberg (FDR) o Bonferroni (FWER) sobre ``pvalues`` y - devuelve, alineado posicion a posicion con la entrada, el p-valor ajustado y + Aplica Benjamini-Hochberg (FDR), Bonferroni (FWER) o Holm-Bonferroni + (FWER, step-down) sobre ``pvalues`` y devuelve, alineado posicion a + posicion con la entrada, el p-valor ajustado y si cada hipotesis se rechaza al nivel ``alpha`` tras la correccion. Las posiciones cuyo valor no sea un p-valor valido (``None``, ``NaN``, fuera de ``[0, 1]`` o no numerico) se conservan en la salida como ``None`` / @@ -53,8 +57,10 @@ def fdr_correction(pvalues: list, alpha: float = 0.05, method: str = "bh") -> di otros valores no validos en posiciones sin test disponible; se propagan como ``None`` en la salida y no cuentan como prueba. alpha: nivel de significancia objetivo tras la correccion (default 0.05). - Para BH es el umbral del FDR; para Bonferroni, del FWER. - method: ``"bh"`` (Benjamini-Hochberg, FDR) o ``"bonferroni"`` (FWER). + Para BH es el umbral del FDR; para Bonferroni y Holm, del FWER. + method: ``"bh"`` (Benjamini-Hochberg, FDR), ``"bonferroni"`` (FWER) o + ``"holm"`` (Holm-Bonferroni, FWER step-down, mas potente que + Bonferroni simple). Returns: dict con las claves: @@ -68,7 +74,7 @@ def fdr_correction(pvalues: list, alpha: float = 0.05, method: str = "bh") -> di n_tests: numero de p-valores validos usados en la correccion (m). n_rejected: numero de hipotesis rechazadas (significativas). alpha: nivel de significancia aplicado (float). - method: metodo aplicado (``"bh"`` o ``"bonferroni"``). + method: metodo aplicado (``"bh"``, ``"bonferroni"`` o ``"holm"``). Casos degenerados (lista vacia, sin p-valores validos o metodo desconocido) anaden ademas una clave ``note`` y devuelven listas @@ -76,7 +82,7 @@ def fdr_correction(pvalues: list, alpha: float = 0.05, method: str = "bh") -> di en las posiciones invalidas). """ method_norm = (method or "").strip().lower() - if method_norm not in {"bh", "bonferroni"}: + if method_norm not in {"bh", "bonferroni", "holm"}: n = len(pvalues) return { "p_values_adjusted": [None] * n, @@ -86,8 +92,8 @@ def fdr_correction(pvalues: list, alpha: float = 0.05, method: str = "bh") -> di "alpha": float(alpha), "method": method, "note": ( - f"metodo desconocido '{method}'; usa 'bh' (Benjamini-Hochberg) " - "o 'bonferroni'" + f"metodo desconocido '{method}'; usa 'bh' (Benjamini-Hochberg), " + "'bonferroni' o 'holm' (Holm-Bonferroni)" ), } @@ -129,6 +135,20 @@ def fdr_correction(pvalues: list, alpha: float = 0.05, method: str = "bh") -> di padj = min(1.0, p * m) adjusted[orig_idx] = padj reject[orig_idx] = padj <= a + elif method_norm == "holm": + # Holm-Bonferroni (step-down). Ordena p ascendente; para el rank k + # (1-indexed) el p ajustado crudo es (m - k + 1) * p_(k). Impon + # monotonicidad acumulada (no decreciente) recorriendo de menor a mayor: + # padj_(k) = max(padj_(k-1), min(1, (m-k+1)*p_(k))), con padj_(0)=0. + order = sorted(valid, key=lambda t: t[1]) # [(orig_idx, p), ...] por p asc + prev = 0.0 + for k in range(1, m + 1): + orig_idx, p = order[k - 1] + raw = min(1.0, (m - k + 1) * p) + padj = max(prev, raw) + prev = padj + adjusted[orig_idx] = padj + reject[orig_idx] = padj <= a else: # Benjamini-Hochberg (step-up). Ordena p ascendente y calcula q-valores # con la monotonicidad acumulada de derecha a izquierda. diff --git a/python/functions/datascience/fdr_correction_test.py b/python/functions/datascience/fdr_correction_test.py index cb6b2bc2..529ed0bd 100644 --- a/python/functions/datascience/fdr_correction_test.py +++ b/python/functions/datascience/fdr_correction_test.py @@ -82,7 +82,8 @@ def test_solo_none_devuelve_note(): def test_metodo_desconocido_devuelve_note(): - out = fdr_correction([0.01, 0.02], method="holm") + # 'holm' ya es un metodo valido (v1.1.0); usamos uno realmente desconocido. + out = fdr_correction([0.01, 0.02], method="sidak") assert "note" in out assert out["n_rejected"] == 0 assert out["reject"] == [False, False] @@ -97,3 +98,66 @@ def test_todos_significativos(): assert bon["n_rejected"] == 3 assert all(bh["reject"]) assert all(bon["reject"]) + + +def test_holm_golden_rechaza_dos_de_cuatro(): + # Holm-Bonferroni (step-down) sobre [0.01, 0.04, 0.03, 0.005], m=4, alpha=0.05. + # Ordenado ascendente: 0.005, 0.01, 0.03, 0.04. + # padj_(1) = 4*0.005 = 0.02 + # padj_(2) = max(0.02, 3*0.01=0.03) = 0.03 + # padj_(3) = max(0.03, 2*0.03=0.06) = 0.06 + # padj_(4) = max(0.06, 1*0.04=0.04) = 0.06 + # Mapeado al orden de entrada [0.01, 0.04, 0.03, 0.005]: + # 0.01 -> 0.03, 0.04 -> 0.06, 0.03 -> 0.06, 0.005 -> 0.02 + out = fdr_correction([0.01, 0.04, 0.03, 0.005], alpha=0.05, method="holm") + assert out["method"] == "holm" + assert out["n_tests"] == 4 + adj = out["p_values_adjusted"] + assert abs(adj[0] - 0.03) < 1e-9 + assert abs(adj[1] - 0.06) < 1e-9 + assert abs(adj[2] - 0.06) < 1e-9 + assert abs(adj[3] - 0.02) < 1e-9 + assert out["reject"] == [True, False, False, True] + assert out["n_rejected"] == 2 + + +def test_holm_entre_bonferroni_y_bh(): + # Holm controla FWER como Bonferroni pero es step-down: rechaza AL MENOS + # tanto como Bonferroni simple, y a lo sumo tanto como BH (FDR, menos + # conservador). Cadena de potencia: bonferroni <= holm <= bh. + pvalues = [0.01, 0.02, 0.04, 0.005] + bon = fdr_correction(pvalues, alpha=0.05, method="bonferroni") + holm = fdr_correction(pvalues, alpha=0.05, method="holm") + bh = fdr_correction(pvalues, alpha=0.05, method="bh") + assert holm["n_rejected"] >= bon["n_rejected"] + assert holm["n_rejected"] <= bh["n_rejected"] + # En este set Holm gana potencia frente a Bonferroni simple (estricto). + assert holm["n_rejected"] > bon["n_rejected"] + + # Un set donde Holm es estrictamente mas conservador que BH. + pvals2 = [0.01, 0.02, 0.03, 0.04] + bon2 = fdr_correction(pvals2, alpha=0.05, method="bonferroni") + holm2 = fdr_correction(pvals2, alpha=0.05, method="holm") + bh2 = fdr_correction(pvals2, alpha=0.05, method="bh") + assert holm2["n_rejected"] >= bon2["n_rejected"] + assert holm2["n_rejected"] < bh2["n_rejected"] + + +def test_none_se_propaga_alineado_holm(): + # None se propaga alineado tambien con holm: la posicion central no cuenta + # como prueba (m=2) y se devuelve como None / False. + out = fdr_correction([0.001, None, 0.9], method="holm") + assert out["n_tests"] == 2 + assert out["p_values_adjusted"][1] is None + assert out["reject"][1] is False + assert out["reject"][0] is True + assert len(out["reject"]) == 3 + + +def test_lista_vacia_holm_devuelve_note(): + out = fdr_correction([], method="holm") + assert out["p_values_adjusted"] == [] + assert out["reject"] == [] + assert out["n_tests"] == 0 + assert out["n_rejected"] == 0 + assert "note" in out diff --git a/python/functions/datascience/preregister_hypothesis.md b/python/functions/datascience/preregister_hypothesis.md new file mode 100644 index 00000000..ea4fc848 --- /dev/null +++ b/python/functions/datascience/preregister_hypothesis.md @@ -0,0 +1,100 @@ +--- +name: preregister_hypothesis +kind: function +lang: py +domain: datascience +version: "1.0.0" +purity: impure +signature: "def preregister_hypothesis(paper_dir: str, hypotheses: dict, analysis_plan: dict) -> dict" +description: "Pre-registra (congela) la hipotesis y el plan de analisis de un paper ANTES de mirar los datos: antidoto al HARKing (Hypothesizing After the Results are Known). Escribe/actualiza /preregistration.md con un frontmatter (paper_slug, frozen_at, content_hash, status) y un cuerpo markdown DETERMINISTA derivado de (hypotheses, analysis_plan) (mismo input -> mismo cuerpo byte a byte, claves ordenadas alfabeticamente). El content_hash es sha256 del cuerpo NORMALIZADO (strip por linea + colapso de blancos), nunca del frontmatter. Una vez status=frozen es INMUTABLE: re-congelar con el mismo contenido es idempotente (no reescribe, devuelve unchanged) y re-congelar con contenido distinto se RECHAZA (no sobrescribe, devuelve error) para que no se pueda ajustar la hipotesis a los resultados. Estilo dict-no-throw: nunca lanza." +tags: [papers, preregistration, reproducibility, anti-harking, python] +params: + - name: paper_dir + desc: "ruta del directorio del paper, p.ej. 'papers/0001-mi-paper'. Debe existir (no se crea aqui). El paper_slug del frontmatter es el basename del dir. Si no existe o no es str -> {status:error, path, note} sin crash ni creacion." + - name: hypotheses + desc: "dict de hipotesis, p.ej. {'h0': 'no hay diferencia ...', 'h1': 'el grupo A > grupo B ...'}. Se renderiza en la seccion '## Hypotheses' con una linea por clave, ordenadas alfabeticamente para determinismo." + - name: analysis_plan + desc: "dict con el plan de analisis, p.ej. {'test': 'welch_t_test', 'effect_size_metric': 'cohens_d', 'decision_rule': 'rechazar H0 si p<0.05 tras Holm y |d|>=0.5', 'planned_n': 100, 'multiple_correction': 'holm'}. Se renderiza en '## Analysis plan' con una linea por clave (ordenadas alfabeticamente). Acepta valores no-str (int, etc.)." +output: "dict dict-no-throw (NUNCA lanza). status='frozen' cuando escribe el archivo por primera vez o congela un draft previo ({status, path, content_hash, frozen_at}). status='unchanged' cuando ya estaba frozen con el mismo content_hash: no reescribe y preserva el archivo byte-identico incl. el frozen_at original ({status, path, content_hash, frozen_at}). status='error' cuando paper_dir no existe, ya esta frozen con un hash distinto (rechazo anti-HARKing, no sobrescribe), inputs invalidos o error de I/O ({status, path, note, [content_hash]})." +uses_functions: [] +uses_types: [] +returns: [] +returns_optional: false +error_type: "error_go_core" +imports: [hashlib] +tested: true +tests: ["test_golden_congela_y_escribe_archivo", "test_idempotente_mismo_input_no_reescribe", "test_inmutabilidad_anti_harking_rechaza_contenido_distinto", "test_error_paper_dir_inexistente_no_crash_no_crea"] +test_file_path: "python/functions/datascience/preregister_hypothesis_test.py" +file_path: "python/functions/datascience/preregister_hypothesis.py" +--- + +## Ejemplo + +```python +import os, tempfile +from datascience import preregister_hypothesis + +# Un directorio de paper que ya existe. +paper_dir = tempfile.mkdtemp(prefix="0001-") + +hypotheses = { + "h0": "no hay diferencia entre el grupo A y el grupo B", + "h1": "el grupo A tiene mayor conversion que el grupo B", +} +analysis_plan = { + "test": "welch_t_test", + "effect_size_metric": "cohens_d", + "decision_rule": "rechazar H0 si p<0.05 tras Holm y |d|>=0.5", + "planned_n": 100, + "multiple_correction": "holm", +} + +# 1) Primera vez: congela y escribe /preregistration.md +r1 = preregister_hypothesis(paper_dir, hypotheses, analysis_plan) +print(r1["status"]) # -> "frozen" +print(r1["content_hash"]) # sha256 del cuerpo + +# 2) Mismo input: idempotente, no reescribe. +r2 = preregister_hypothesis(paper_dir, hypotheses, analysis_plan) +print(r2["status"]) # -> "unchanged" + +# 3) Cambiar la hipotesis tras congelar (HARKing): rechazado, archivo intacto. +r3 = preregister_hypothesis(paper_dir, {"h0": "...", "h1": "otra cosa"}, analysis_plan) +print(r3["status"]) # -> "error" +``` + +## Cuando usarla + +Llamala al ARRANCAR el analisis de un paper, antes de tocar los datos, para +dejar por escrito (y firmado por hash) que vas a probar y como vas a decidir. +Es el primer paso de un flujo reproducible: pre-registras la hipotesis y el plan +(`test`, `effect_size_metric`, `decision_rule`, `planned_n`, +`multiple_correction`), y solo despues corres el analisis y comparas con lo +pre-registrado. Si mas tarde el analisis "descubre" otra hipotesis que encaja +mejor con los datos, el pre-registro congelado deja en evidencia el cambio: no se +puede reescribir. Combinala con `effect_size_cohens_d` y `fdr_correction` para +cerrar el plan declarado (effect size + correccion de multiples comparaciones). + +## Gotchas + +- **Inmutabilidad (el corazon)**: una vez `status: frozen`, el pre-registro NO se + puede editar. Re-congelar con el MISMO contenido es idempotente (`unchanged`, + no reescribe, preserva incluso el `frozen_at` original). Re-congelar con + contenido DISTINTO devuelve `error` y deja el archivo intacto: asi se mata el + HARKing. Para cambiar de verdad la hipotesis hay que borrar el archivo a mano y + asumir explicitamente que ya no es un pre-registro valido. +- **dict-no-throw**: la funcion NUNCA lanza. Cualquier error previsible + (directorio inexistente, inputs no-dict, fallo de I/O, excepcion inesperada) se + captura y se devuelve como `{"status": "error", "note": ...}`. Siempre incluye + `path` (la ruta esperada del `preregistration.md`). +- **El hash es SOLO del cuerpo, nunca del frontmatter**: el frontmatter contiene + el propio `content_hash` y el `frozen_at` (timestamp), asi que incluirlos en el + hash seria circular y romperia la idempotencia. El cuerpo se normaliza antes de + hashear (strip por linea + colapso de lineas en blanco + strip final): cambios + irrelevantes de whitespace no alteran el hash, pero cambios de contenido SI. +- **Determinismo**: el cuerpo se genera con las claves de `hypotheses` y + `analysis_plan` ordenadas alfabeticamente, de modo que el orden de insercion del + dict no afecta al hash. Mismo `(hypotheses, analysis_plan)` -> mismo cuerpo y + mismo hash, byte a byte. +- **No crea el directorio del paper**: si `paper_dir` no existe, devuelve `error` + sin crear nada (ni el dir ni el archivo). diff --git a/python/functions/datascience/preregister_hypothesis.py b/python/functions/datascience/preregister_hypothesis.py new file mode 100644 index 00000000..0be8cc2b --- /dev/null +++ b/python/functions/datascience/preregister_hypothesis.py @@ -0,0 +1,202 @@ +"""Congela (pre-registra) la hipotesis y el plan de analisis de un paper. + +Anti-HARKing (Hypothesizing After the Results are Known): el pre-registro fija +la hipotesis y el plan de analisis ANTES de mirar los datos. Una vez congelado +(``status: frozen``) es INMUTABLE: cualquier intento posterior de re-congelar con +un contenido distinto se RECHAZA en vez de sobrescribir, de modo que no se puede +"ajustar" la hipotesis a los resultados despues de verlos. + +Escribe/actualiza ``/preregistration.md`` con un frontmatter +(``paper_slug``, ``frozen_at``, ``content_hash``, ``status``) y un cuerpo +markdown DETERMINISTA derivado de ``(hypotheses, analysis_plan)``. + +Estilo dict-no-throw: NUNCA lanza; cualquier error previsible se captura y se +devuelve como ``{"status": "error", "note": ...}``. +""" + +import hashlib +import os +from datetime import datetime, timezone + + +def _build_body(hypotheses: dict, analysis_plan: dict) -> str: + """Construye el cuerpo markdown del pre-registro de forma DETERMINISTA. + + Mismo ``(hypotheses, analysis_plan)`` -> mismo cuerpo byte a byte. Las claves + se ordenan alfabeticamente para no depender del orden de insercion del dict. + """ + lines = ["## Hypotheses", ""] + for k in sorted(hypotheses.keys()): + lines.append(f"- **{k}**: {hypotheses[k]}") + lines.append("") + lines.append("## Analysis plan") + lines.append("") + for k in sorted(analysis_plan.keys()): + lines.append(f"- **{k}**: {analysis_plan[k]}") + return "\n".join(lines) + + +def _normalize(body: str) -> str: + """Normaliza el cuerpo para el hash: strip por linea + colapsa blancos. + + Cambios irrelevantes de whitespace (espacios al final, dobles lineas en + blanco) no alteran el hash; cambios de contenido SI. Esto hace el hash + robusto sin perder la capacidad de detectar ediciones reales. + """ + out = [] + prev_blank = False + for raw in body.splitlines(): + line = raw.strip() + if line == "": + if prev_blank: + continue + prev_blank = True + else: + prev_blank = False + out.append(line) + return "\n".join(out).strip() + + +def _content_hash(body: str) -> str: + """sha256 hex del cuerpo NORMALIZADO (nunca del frontmatter).""" + return hashlib.sha256(_normalize(body).encode("utf-8")).hexdigest() + + +def _parse_frontmatter(text: str) -> dict: + """Parsea el frontmatter ``--- ... ---`` simple (key: value) de un .md.""" + if not text.startswith("---"): + return {} + parts = text.split("---", 2) + if len(parts) < 3: + return {} + fm = {} + for line in parts[1].splitlines(): + line = line.strip() + if not line or ":" not in line: + continue + key, _, value = line.partition(":") + fm[key.strip()] = value.strip() + return fm + + +def _render_file(slug: str, frozen_at: str, content_hash: str, body: str) -> str: + """Compone el archivo completo: frontmatter frozen + cuerpo.""" + return ( + "---\n" + f"paper_slug: {slug}\n" + f"frozen_at: {frozen_at}\n" + f"content_hash: {content_hash}\n" + "status: frozen\n" + "---\n" + "\n" + f"{body}\n" + ) + + +def preregister_hypothesis(paper_dir: str, hypotheses: dict, analysis_plan: dict) -> dict: + """Congela la hipotesis y el plan de analisis de un paper (anti-HARKing). + + Escribe ``/preregistration.md`` con frontmatter ``status: frozen`` + y un cuerpo markdown determinista. Una vez congelado es inmutable. + + Args: + paper_dir: ruta del directorio del paper (p.ej. ``"papers/0001-mi-paper"``). + El ``paper_slug`` es el basename del directorio. Debe existir. + hypotheses: dict de hipotesis, p.ej. + ``{"h0": "no hay diferencia ...", "h1": "grupo A > grupo B ..."}``. + analysis_plan: dict con el plan, p.ej. + ``{"test": "welch_t_test", "effect_size_metric": "cohens_d", + "decision_rule": "...", "planned_n": 100, "multiple_correction": "holm"}``. + + Returns: + dict dict-no-throw (NUNCA lanza). Claves segun el caso: + - frozen: {"status": "frozen", "path", "content_hash", "frozen_at"} + - unchanged: {"status": "unchanged", "path", "content_hash", "frozen_at"} + - error: {"status": "error", "path", "note", ...} + """ + expected_path = os.path.join(paper_dir, "preregistration.md") + try: + # 1) El directorio del paper debe existir; no se crea aqui. + if not isinstance(paper_dir, str) or not os.path.isdir(paper_dir): + return { + "status": "error", + "path": expected_path, + "note": f"paper_dir no existe: {paper_dir}", + } + + if not isinstance(hypotheses, dict) or not isinstance(analysis_plan, dict): + return { + "status": "error", + "path": expected_path, + "note": "hypotheses y analysis_plan deben ser dict", + } + + slug = os.path.basename(os.path.normpath(paper_dir)) + + # 2) + 3) Cuerpo determinista y su hash (solo del cuerpo, no del frontmatter). + body = _build_body(hypotheses, analysis_plan) + new_hash = _content_hash(body) + + # 5) Logica de escritura. + if os.path.exists(expected_path): + existing = "" + try: + with open(expected_path, "r", encoding="utf-8") as fh: + existing = fh.read() + except OSError as exc: + return { + "status": "error", + "path": expected_path, + "note": f"no se pudo leer el pre-registro existente: {exc}", + } + fm = _parse_frontmatter(existing) + old_status = fm.get("status", "") + old_hash = fm.get("content_hash", "") + old_frozen_at = fm.get("frozen_at", "") + + if old_status == "frozen": + if old_hash == new_hash: + # Idempotente: mismo contenido ya congelado. No se reescribe. + return { + "status": "unchanged", + "path": expected_path, + "content_hash": new_hash, + "frozen_at": old_frozen_at, + } + # Inmutabilidad: ya congelado con OTRO hash -> se rechaza (anti-HARKing). + return { + "status": "error", + "path": expected_path, + "content_hash": new_hash, + "note": ( + "pre-registro inmutable: ya esta congelado (frozen) con un " + "hash distinto; un pre-registro no se puede editar tras " + "congelarse" + ), + } + # status != "frozen" (p.ej. draft) -> se congela ahora. + + # Archivo nuevo o draft existente: congelar con timestamp actual. + frozen_at = datetime.now(timezone.utc).strftime("%Y-%m-%dT%H:%M:%SZ") + file_text = _render_file(slug, frozen_at, new_hash, body) + try: + with open(expected_path, "w", encoding="utf-8") as fh: + fh.write(file_text) + except OSError as exc: + return { + "status": "error", + "path": expected_path, + "note": f"no se pudo escribir el pre-registro: {exc}", + } + return { + "status": "frozen", + "path": expected_path, + "content_hash": new_hash, + "frozen_at": frozen_at, + } + except Exception as exc: # noqa: BLE001 - dict-no-throw: nunca propagar. + return { + "status": "error", + "path": expected_path, + "note": f"error inesperado: {exc}", + } diff --git a/python/functions/datascience/preregister_hypothesis_test.py b/python/functions/datascience/preregister_hypothesis_test.py new file mode 100644 index 00000000..93c71e3e --- /dev/null +++ b/python/functions/datascience/preregister_hypothesis_test.py @@ -0,0 +1,99 @@ +"""Tests para preregister_hypothesis (pre-registro inmutable, anti-HARKing). + +Importa el modulo hoja directamente (`preregister_hypothesis`) para no depender +de que el paquete reexporte la funcion en su __init__ (lo integra el orquestador +al cerrar el grupo papers). El pytest del repo resuelve el modulo hoja por su +nombre directo. + +Todos los tests son hermeticos y deterministas: usan el fixture `tmp_path` de +pytest; NUNCA escriben en `papers/`. +""" + +from preregister_hypothesis import preregister_hypothesis + + +def _parse_frontmatter(text: str) -> dict: + parts = text.split("---", 2) + fm = {} + for line in parts[1].splitlines(): + line = line.strip() + if not line or ":" not in line: + continue + key, _, value = line.partition(":") + fm[key.strip()] = value.strip() + return fm + + +HYP = {"h0": "no hay diferencia entre A y B", "h1": "el grupo A > grupo B"} +PLAN = { + "test": "welch_t_test", + "effect_size_metric": "cohens_d", + "decision_rule": "rechazar H0 si p<0.05 tras Holm y |d|>=0.5", + "planned_n": 100, + "multiple_correction": "holm", +} + + +def test_golden_congela_y_escribe_archivo(tmp_path): + paper = tmp_path / "0001-x" + paper.mkdir() + + res = preregister_hypothesis(str(paper), HYP, PLAN) + + assert res["status"] == "frozen" + pre = paper / "preregistration.md" + assert pre.exists() + + text = pre.read_text(encoding="utf-8") + fm = _parse_frontmatter(text) + assert fm["status"] == "frozen" + assert fm["paper_slug"] == "0001-x" + assert fm["content_hash"] # no vacio + assert fm["frozen_at"] # no vacio + assert res["content_hash"] == fm["content_hash"] + assert res["frozen_at"] == fm["frozen_at"] + + +def test_idempotente_mismo_input_no_reescribe(tmp_path): + paper = tmp_path / "0001-x" + paper.mkdir() + pre = paper / "preregistration.md" + + first = preregister_hypothesis(str(paper), HYP, PLAN) + assert first["status"] == "frozen" + bytes_before = pre.read_bytes() + + second = preregister_hypothesis(str(paper), HYP, PLAN) + assert second["status"] == "unchanged" + # Mismo hash y frozen_at original preservado. + assert second["content_hash"] == first["content_hash"] + assert second["frozen_at"] == first["frozen_at"] + # El archivo NO cambio byte a byte (incl. frozen_at). + assert pre.read_bytes() == bytes_before + + +def test_inmutabilidad_anti_harking_rechaza_contenido_distinto(tmp_path): + paper = tmp_path / "0001-x" + paper.mkdir() + pre = paper / "preregistration.md" + + preregister_hypothesis(str(paper), HYP, PLAN) + bytes_frozen = pre.read_bytes() + + # Intento de re-congelar con una hipotesis DISTINTA (HARKing) -> rechazado. + hyp_tramposo = {"h0": "no hay diferencia", "h1": "el grupo B > grupo A (cambiado tras ver datos)"} + res = preregister_hypothesis(str(paper), hyp_tramposo, PLAN) + + assert res["status"] == "error" + # Asercion mas importante: el archivo en disco SIGUE siendo el original. + assert pre.read_bytes() == bytes_frozen + + +def test_error_paper_dir_inexistente_no_crash_no_crea(tmp_path): + missing = tmp_path / "no-existe" + res = preregister_hypothesis(str(missing), HYP, PLAN) + + assert res["status"] == "error" + # No se creo el directorio ni el archivo. + assert not missing.exists() + assert not (missing / "preregistration.md").exists()