 codex/nail-kernel-centroid-ambiguity.py

"""Two original normalized material-age kernels with the same centroid.

Python 3.14, NumPy 2.3.5, Matplotlib 3.10.7. Writes PNG into cwd.
The plotted kernels are hypothetical reference-crossing ages, not exposure data.
"""
from pathlib import Path

import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle


def moments(windows):
    """Exact mean and variance for (lower age, upper age, mass fraction)."""
    assert abs(sum(weight for _, _, weight in windows) - 1.0) < 1e-12
    assert all(upper > lower and weight >= 0 for lower, upper, weight in windows)
    mean = sum(weight * (lower + upper) / 2 for lower, upper, weight in windows)
    second = sum(weight * (lower**2 + lower * upper + upper**2) / 3
                 for lower, upper, weight in windows)
    return mean, second - mean**2


def main():
    cases = [([(0, 20, 1.0)], "Uniform material ages", "#3979a8", 0.08),
             ([(0, 2, 0.5), (18, 20, 0.5)], "Two pooled, separated segments", "#438354", 0.32)]
    fig, axes = plt.subplots(2, 1, figsize=(7.2, 4.8), dpi=100, sharex=True,
                             facecolor="white")
    for ax, (windows, label, color, ymax) in zip(axes, cases):
        ax.set_facecolor("white")
        mean, variance = moments(windows)
        for lower, upper, weight in windows:
            ax.add_patch(Rectangle((lower, 0), upper - lower,
                                   weight / (upper - lower), facecolor=color,
                                   edgecolor=color, alpha=0.65))
        ax.axvline(mean, color="#9a4a35", ls="--", lw=1.7)
        ax.set(xlim=(-0.5, 20.5), ylim=(0, ymax), ylabel="Mass fraction / day")
        ax.set_title(f"{label} — centroid 10 days; variance {variance:.2f} days²", fontsize=10)
        ax.grid(axis="y", alpha=0.18)
    axes[1].text(10, 0.17, "No sampled material\nat the centroid age", ha="center", fontsize=9,
                 bbox={"facecolor": "white", "edgecolor": "none", "alpha": 0.95})
    axes[1].set_xlabel("Material reference-crossing age before collection (days)")
    axes[1].set_xticks([0, 2, 5, 10, 15, 18, 20])
    fig.suptitle("The same centroid can summarize very different material histories", fontsize=11)
    fig.text(0.5, 0.015, "Original normalized kernels: each shaded distribution has total mass fraction 1.\n"
             "These are material-age models, not measured exposure histories or clinical dating rules.",
             ha="center", fontsize=8)
    fig.tight_layout(rect=(0, 0.07, 1, 0.96))
    fig.savefig(Path("nail-kernel-centroid-ambiguity.png"), facecolor="white")
    plt.close(fig)


if __name__ == "__main__":
    main()