From 22787da898003458f0d394a22e5591bb7c6c1986 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 21 Aug 2024 15:53:39 +0000 Subject: [PATCH] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- doc/notebooks/weighted_histograms.ipynb | 115 +++++++++++++++++++----- 1 file changed, 94 insertions(+), 21 deletions(-) diff --git a/doc/notebooks/weighted_histograms.ipynb b/doc/notebooks/weighted_histograms.ipynb index 4be26d73..8b22fcce 100644 --- a/doc/notebooks/weighted_histograms.ipynb +++ b/doc/notebooks/weighted_histograms.ipynb @@ -109,6 +109,7 @@ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", + "\n", "def log_or_zero(x):\n", " x = np.atleast_1d(x).copy()\n", " ma = x > 0\n", @@ -116,18 +117,22 @@ " x[~ma] = 0\n", " return x\n", "\n", + "\n", "def l1(w, w2, mu):\n", " s = np.abs(w) / w2\n", " return 2 * s * (w * (log_or_zero(w) - log_or_zero(mu)) + mu - w)\n", "\n", + "\n", "def l2(w, w2, mu):\n", " s = np.abs(w) / w2\n", " return 2 * s * (w * log_or_zero(np.abs(w)) - w * log_or_zero(mu) + mu - w)\n", "\n", + "\n", "def l3(w, w2, mu):\n", " s = w / w2\n", " return 2 * s * (w * (log_or_zero(np.abs(w)) - log_or_zero(mu)) + mu - w)\n", "\n", + "\n", "w2 = 1\n", "fig, ax = plt.subplots(1, 3, figsize=(10, 4))\n", "for i, (axi, w) in enumerate(zip(ax, (-10, -5, 5))):\n", @@ -190,13 +195,23 @@ "x = expon.rvs(size=npoints, random_state=rng)\n", "w = norm.rvs(0.1, x, size=len(x), random_state=rng)\n", "\n", - "h = bh.Histogram(bh.axis.Regular(bins, np.min(x), np.max(x)), storage=bh.storage.Weight())\n", + "h = bh.Histogram(\n", + " bh.axis.Regular(bins, np.min(x), np.max(x)), storage=bh.storage.Weight()\n", + ")\n", "h.fill(x, weight=w)\n", "\n", "plt.stairs(h.values(), h.axes[0].edges)\n", "ma = h.values() > 0\n", - "plt.errorbar(h.axes[0].centers[ma], h.values()[ma], h.variances()[ma] ** 0.5, fmt=\"o\", color=\"C0\")\n", - "plt.errorbar(h.axes[0].centers[~ma], h.values()[~ma], h.variances()[~ma] ** 0.5, fmt=\"s\", color=\"C3\")\n", + "plt.errorbar(\n", + " h.axes[0].centers[ma], h.values()[ma], h.variances()[ma] ** 0.5, fmt=\"o\", color=\"C0\"\n", + ")\n", + "plt.errorbar(\n", + " h.axes[0].centers[~ma],\n", + " h.values()[~ma],\n", + " h.variances()[~ma] ** 0.5,\n", + " fmt=\"s\",\n", + " color=\"C3\",\n", + ")\n", "plt.axhline(0, ls=\":\", color=\"0.5\", zorder=0)\n", "xm = np.linspace(0, h.axes[0].edges[-1])\n", "plt.plot(xm, expon.pdf(xm) * npoints * np.mean(w) * h.axes[0].widths[0]);" @@ -931,6 +946,7 @@ "def model1(x, n, lambd):\n", " return n * expon(0, lambd).cdf(x)\n", "\n", + "\n", "c1 = ExtendedBinnedNLL(np.transpose((n, vn)), xe, model1)\n", "m = Minuit(c1, sum(n), 1)\n", "m.limits = (0, None)\n", @@ -1621,6 +1637,7 @@ "def model2(x, lambd):\n", " return expon(0, lambd).cdf(x)\n", "\n", + "\n", "c2 = BinnedNLL(np.transpose((n, vn)), xe, model2)\n", "m = Minuit(c2, 1)\n", "m.limits = (0, None)\n", @@ -1648,7 +1665,9 @@ " x = expon.rvs(size=rng.poisson(npoints), random_state=rng)\n", " w = norm.rvs(0.1, x, size=len(x), random_state=rng)\n", "\n", - " h = bh.Histogram(bh.axis.Regular(bins, np.min(x), np.max(x)), storage=bh.storage.Weight())\n", + " h = bh.Histogram(\n", + " bh.axis.Regular(bins, np.min(x), np.max(x)), storage=bh.storage.Weight()\n", + " )\n", " h.fill(x, weight=w)\n", " xe = h.axes[0].edges\n", " n = h.values()\n", @@ -1665,7 +1684,17 @@ " m2.limits = (0, None)\n", " m2.migrad()\n", "\n", - " return ntot, m1.valid, m1.values[0], m1.values[1], m1.fval, m2.valid, m2.values[0], m2.fval\n", + " return (\n", + " ntot,\n", + " m1.valid,\n", + " m1.values[0],\n", + " m1.values[1],\n", + " m1.fval,\n", + " m2.valid,\n", + " m2.values[0],\n", + " m2.fval,\n", + " )\n", + "\n", "\n", "result = Parallel(n_jobs=8)(delayed(run)(seed) for seed in range(1000))\n", "ntot, valid1, ntot1, lambd1, minval1, valid2, lambd2, minval2 = np.transpose(result)" @@ -1700,13 +1729,25 @@ "source": [ "plt.figure()\n", "plt.title(\"$\\\\lambda$ estimate\")\n", - "plt.hist(lambd1, bins=20, alpha=0.5, label=f\"ExtendedBinnedNLL\\nvalid={np.mean(valid1)*100:.0f}%\\nmean = {np.mean(lambd1):.2f}\\nstd.dev. = {np.std(lambd1):.2f}\")\n", - "plt.hist(lambd2, bins=20, alpha=0.5, label=f\"BinnedNLL\\nvalid={np.mean(valid2)*100:.0f}%\\nmean = {np.mean(lambd2):.2f}\\nstd.dev. = {np.std(lambd2):.2f}\")\n", + "plt.hist(\n", + " lambd1,\n", + " bins=20,\n", + " alpha=0.5,\n", + " label=f\"ExtendedBinnedNLL\\nvalid={np.mean(valid1)*100:.0f}%\\nmean = {np.mean(lambd1):.2f}\\nstd.dev. = {np.std(lambd1):.2f}\",\n", + ")\n", + "plt.hist(\n", + " lambd2,\n", + " bins=20,\n", + " alpha=0.5,\n", + " label=f\"BinnedNLL\\nvalid={np.mean(valid2)*100:.0f}%\\nmean = {np.mean(lambd2):.2f}\\nstd.dev. = {np.std(lambd2):.2f}\",\n", + ")\n", "plt.legend()\n", "\n", "plt.figure()\n", "plt.hist(ntot1 / ntot - 1)\n", - "plt.title(f\"ExtendedBinnedNLL: $n_\\\\mathrm{{tot}}$ estimate\\n(mean-truth)/truth = {np.mean(ntot1) / np.mean(ntot) - 1:.3f}\")\n", + "plt.title(\n", + " f\"ExtendedBinnedNLL: $n_\\\\mathrm{{tot}}$ estimate\\n(mean-truth)/truth = {np.mean(ntot1) / np.mean(ntot) - 1:.3f}\"\n", + ")\n", "plt.xlabel(\"(estimate - truth) / truth\")\n", "plt.xlim(-0.2, 0.2);" ] @@ -1748,14 +1789,18 @@ "plt.sca(ax[0])\n", "plt.hist(minval1, bins=50, density=True)\n", "x = np.linspace(0, 100)\n", - "plt.plot(x, chi2(bins-2).pdf(x))\n", - "plt.title(f\"ExtendedBinnedNLL minimum value\\nndf = {bins-2} mean = {np.mean(minval1):.2f} median = {np.median(minval1):.2f}\")\n", + "plt.plot(x, chi2(bins - 2).pdf(x))\n", + "plt.title(\n", + " f\"ExtendedBinnedNLL minimum value\\nndf = {bins-2} mean = {np.mean(minval1):.2f} median = {np.median(minval1):.2f}\"\n", + ")\n", "\n", "plt.sca(ax[1])\n", "plt.hist(minval2, bins=50, density=True)\n", "x = np.linspace(0, 100)\n", - "plt.plot(x, chi2(bins-2).pdf(x))\n", - "plt.title(f\"BinnedNLL minimum value\\nndf = {bins-2} mean = {np.mean(minval2):.2f} median = {np.median(minval2):.2f}\");" + "plt.plot(x, chi2(bins - 2).pdf(x))\n", + "plt.title(\n", + " f\"BinnedNLL minimum value\\nndf = {bins-2} mean = {np.mean(minval2):.2f} median = {np.median(minval2):.2f}\"\n", + ");" ] }, { @@ -1784,7 +1829,9 @@ " x = expon.rvs(size=rng.poisson(npoints), random_state=rng)\n", " w = rng.normal(0.1, 0.01, size=len(x))\n", "\n", - " h = bh.Histogram(bh.axis.Regular(bins, np.min(x), np.max(x)), storage=bh.storage.Weight())\n", + " h = bh.Histogram(\n", + " bh.axis.Regular(bins, np.min(x), np.max(x)), storage=bh.storage.Weight()\n", + " )\n", " h.fill(x, weight=w)\n", " xe = h.axes[0].edges\n", " n = h.values()\n", @@ -1801,7 +1848,17 @@ " m2.limits = (0, None)\n", " m2.migrad()\n", "\n", - " return ntot, m1.valid, m1.values[0], m1.values[1], m1.fval, m2.valid, m2.values[0], m2.fval\n", + " return (\n", + " ntot,\n", + " m1.valid,\n", + " m1.values[0],\n", + " m1.values[1],\n", + " m1.fval,\n", + " m2.valid,\n", + " m2.values[0],\n", + " m2.fval,\n", + " )\n", + "\n", "\n", "result = Parallel(n_jobs=8)(delayed(run)(seed) for seed in range(1000))\n", "ntot, valid1, ntot1, lambd1, minval1, valid2, lambd2, minval2 = np.transpose(result)" @@ -1836,13 +1893,25 @@ "source": [ "plt.figure()\n", "plt.title(\"$\\\\lambda$ estimate\")\n", - "plt.hist(lambd1, bins=20, alpha=0.5, label=f\"ExtendedBinnedNLL\\nvalid={np.mean(valid1)*100:.0f}%\\nmean = {np.mean(lambd1):.2f}\\nstd.dev. = {np.std(lambd1):.3f}\")\n", - "plt.hist(lambd2, bins=20, alpha=0.5, label=f\"BinnedNLL\\nvalid={np.mean(valid2)*100:.0f}%\\nmean = {np.mean(lambd2):.2f}\\nstd.dev. = {np.std(lambd2):.3f}\")\n", + "plt.hist(\n", + " lambd1,\n", + " bins=20,\n", + " alpha=0.5,\n", + " label=f\"ExtendedBinnedNLL\\nvalid={np.mean(valid1)*100:.0f}%\\nmean = {np.mean(lambd1):.2f}\\nstd.dev. = {np.std(lambd1):.3f}\",\n", + ")\n", + "plt.hist(\n", + " lambd2,\n", + " bins=20,\n", + " alpha=0.5,\n", + " label=f\"BinnedNLL\\nvalid={np.mean(valid2)*100:.0f}%\\nmean = {np.mean(lambd2):.2f}\\nstd.dev. = {np.std(lambd2):.3f}\",\n", + ")\n", "plt.legend()\n", "\n", "plt.figure()\n", "plt.hist(ntot1 / ntot - 1, bins=20)\n", - "plt.title(f\"ExtendedBinnedNLL: $n_\\\\mathrm{{tot}}$ estimate\\n(mean-truth)/truth = {np.mean(ntot1) / np.mean(ntot) - 1:.3f}\")\n", + "plt.title(\n", + " f\"ExtendedBinnedNLL: $n_\\\\mathrm{{tot}}$ estimate\\n(mean-truth)/truth = {np.mean(ntot1) / np.mean(ntot) - 1:.3f}\"\n", + ")\n", "plt.xlabel(\"(estimate - truth) / truth\");" ] }, @@ -1879,14 +1948,18 @@ "plt.sca(ax[0])\n", "plt.hist(minval1, bins=50, density=True)\n", "x = np.linspace(0, 100)\n", - "plt.plot(x, chi2(bins-2).pdf(x))\n", - "plt.title(f\"ExtendedBinnedNLL minimum value\\nndf = {bins-2} mean = {np.mean(minval1):.2f} median = {np.median(minval1):.2f}\")\n", + "plt.plot(x, chi2(bins - 2).pdf(x))\n", + "plt.title(\n", + " f\"ExtendedBinnedNLL minimum value\\nndf = {bins-2} mean = {np.mean(minval1):.2f} median = {np.median(minval1):.2f}\"\n", + ")\n", "\n", "plt.sca(ax[1])\n", "plt.hist(minval2, bins=50, density=True)\n", "x = np.linspace(0, 100)\n", - "plt.plot(x, chi2(bins-2).pdf(x))\n", - "plt.title(f\"BinnedNLL minimum value\\nndf = {bins-2} mean = {np.mean(minval2):.2f} median = {np.median(minval2):.2f}\");" + "plt.plot(x, chi2(bins - 2).pdf(x))\n", + "plt.title(\n", + " f\"BinnedNLL minimum value\\nndf = {bins-2} mean = {np.mean(minval2):.2f} median = {np.median(minval2):.2f}\"\n", + ");" ] }, {