Use lower and upper to match scipy terms

Split warning for readability Add test for rtol warning Remove tmpdir fixture as not needed for these tests given no writing of output Add check for return_results More verbose fix: Correct concatenate lists instead of adding float to all list elements Test bounds expansion test: Update test_plot_results_no_axis baseline image (#2009) * matplotlib v3.6.0 results in a slightly different baseline image than matplotlib v3.5.x, so regenerate the baseline image using matplotlib v3.6.0 with `pytest --mpl-generate-path=tests/contrib/baseline tests/contrib/test_viz.py`. * Mark the test_plot_results_no_axis test as xfail for Python 3.7 as matplotlib v3.6.0 is Python 3.8+ and so the image is guaranteed to be different as Python 3.7 runtimes will install matplotlib v3.5.x. Add upperlimit_fixed_scan to API docs Add return_results test move to test_upperlimit_with_kwargs Move the pop out before evaluation to make everything very clean and clear Note what scan Rename to auto_scan docs: fix link Provide better coverage and use np.allclose docs: Add Beojan Stanislaus to contributor list change auto_scan to toms748_scan rename fixed_scan to linear_grid_scan Make intervals module and change API to upper_limit Rename to pyhf.infer.intervals.upper_limits get upper_limits.upper_limit working Also bring along old API limit to just upper_limit by default Rearrange feat: Add internal API to warn of deprecation and future removal * Add internal API pyhf.exceptions._deprecated_api_warning to alert users to API deprecation by raising a subclass of DeprecationWarning and future removal. * Add test for pyhf.exceptions._deprecated_api_warning to ensure it gets picked up as DeprecationWarning. Note deprecated API Seperate into condifence intervals section fix: Use function scope to avoid altering hypotest_args fixture Make test name explicit Use deprecated Sphinx note Add versionadded directives feat: Add internal API to warn of deprecation and future removal (#2012) * Add internal API pyhf.exceptions._deprecated_api_warning to alert users to API deprecation by raising a subclass of DeprecationWarning and future removal. * Add test for pyhf.exceptions._deprecated_api_warning to ensure it gets picked up as DeprecationWarning. Update lower bound on scipy as toms748 added in scipy v1.2.0
scikit-hep · Sep 20, 2022 · dffb458 · dffb458
1 parent 19a2ee1
commit dffb458
Show file tree

Hide file tree

Showing 15 changed files with 258 additions and 76 deletions.
diff --git a/docs/api.rst b/docs/api.rst
@@ -157,10 +157,19 @@ Fits and Tests
    mle.fit
    mle.fixed_poi_fit
    hypotest
-   intervals.upperlimit
-   intervals.upperlimit_auto
    utils.all_pois_floating
 
+Confidence Intervals
+~~~~~~~~~~~~~~~~~~~~
+
+.. autosummary::
+   :toctree: _generated/
+   :nosignatures:
+
+   intervals.upper_limits.upper_limit
+   intervals.upper_limits.toms748_scan
+   intervals.upper_limits.linear_grid_scan
+   intervals.upperlimit
 
 Schema
 ------

diff --git a/docs/contributors.rst b/docs/contributors.rst
@@ -29,3 +29,4 @@ Contributors include:
 - Aryan Roy
 - Jerry Ling
 - Nathan Simpson
+- Beojan Stanislaus
diff --git a/docs/examples/notebooks/ShapeFactor.ipynb b/docs/examples/notebooks/ShapeFactor.ipynb
@@ -176,7 +176,11 @@
     }
    ],
    "source": [
-    "obs_limit, exp_limits, (poi_tests, tests) = pyhf.infer.intervals.upperlimit(\n",
+    "(\n",
+    "    obs_limit,\n",
+    "    exp_limits,\n",
+    "    (poi_tests, tests),\n",
+    ") = pyhf.infer.intervals.upper_limits.upper_limit(\n",
     "    data, pdf, np.linspace(0, 5, 61), level=0.05, return_results=True\n",
     ")"
    ]

diff --git a/docs/examples/notebooks/binderexample/StatisticalAnalysis.ipynb b/docs/examples/notebooks/binderexample/StatisticalAnalysis.ipynb
@@ -2295,7 +2295,11 @@
    "outputs": [],
    "source": [
     "mu_tests = np.linspace(0, 1, 16)\n",
-    "obs_limit, exp_limits, (poi_tests, tests) = pyhf.infer.intervals.upperlimit(\n",
+    "(\n",
+    "    obs_limit,\n",
+    "    exp_limits,\n",
+    "    (poi_tests, tests),\n",
+    ") = pyhf.infer.intervals.upper_limits.upper_limit(\n",
     "    data, pdf, mu_tests, level=0.05, return_results=True\n",
     ")"
    ]

diff --git a/docs/examples/notebooks/multiBinPois.ipynb b/docs/examples/notebooks/multiBinPois.ipynb
@@ -94,7 +94,11 @@
     "init_pars = model.config.suggested_init()\n",
     "par_bounds = model.config.suggested_bounds()\n",
     "\n",
-    "obs_limit, exp_limits, (poi_tests, tests) = pyhf.infer.intervals.upperlimit(\n",
+    "(\n",
+    "    obs_limit,\n",
+    "    exp_limits,\n",
+    "    (poi_tests, tests),\n",
+    ") = pyhf.infer.intervals.upper_limits.upper_limit(\n",
     "    data, model, np.linspace(0, 5, 61), level=0.05, return_results=True\n",
     ")"
    ]

diff --git a/setup.cfg b/setup.cfg
@@ -35,7 +35,7 @@ package_dir =
 include_package_data = True
 python_requires = >=3.7
 install_requires =
-    scipy>=1.1.0 # requires numpy, which is required by pyhf and tensorflow
+    scipy>=1.2.0 # requires numpy, which is required by pyhf and tensorflow
     click>=8.0.0  # for console scripts
     tqdm>=4.56.0  # for readxml
     jsonschema>=4.15.0  # for utils

diff --git a/src/pyhf/exceptions/__init__.py b/src/pyhf/exceptions/__init__.py
@@ -1,4 +1,5 @@
 import sys
+from warnings import warn
 
 __all__ = [
     "FailedMinimization",
@@ -175,3 +176,15 @@ def __init__(self, result):
             result, 'message', "Unknown failure. See fit result for more details."
         )
         super().__init__(message)
+
+
+# Deprecated APIs
+def _deprecated_api_warning(
+    deprecated_api, new_api, deprecated_release, remove_release
+):
+    warn(
+        f"{deprecated_api} is deprecated in favor of {new_api} as of pyhf v{deprecated_release} and will be removed in pyhf v{remove_release}."
+        + f" Please use {new_api}.",
+        DeprecationWarning,
+        stacklevel=3,  # Raise to user level
+    )
diff --git a/src/pyhf/infer/intervals/__init__.py b/src/pyhf/infer/intervals/__init__.py
@@ -0,0 +1,30 @@
+"""Interval estimation"""
+import pyhf.infer.intervals.upper_limits
+
+__all__ = ["upper_limits.upper_limit"]
+
+
+def __dir__():
+    return __all__
+
+
+def upperlimit(
+    data, model, scan=None, level=0.05, return_results=False, **hypotest_kwargs
+):
+    """
+    .. deprecated:: 0.7.0
+       Use :func:`~pyhf.infer.intervals.upper_limits.upper_limit` instead.
+    .. warning:: :func:`~pyhf.infer.intervals.upperlimit` will be removed in
+     ``pyhf`` ``v0.9.0``.
+    """
+    from pyhf.exceptions import _deprecated_api_warning
+
+    _deprecated_api_warning(
+        "pyhf.infer.intervals.upperlimit",
+        "pyhf.infer.intervals.upper_limits.upper_limit",
+        "0.7.0",
+        "0.9.0",
+    )
+    return pyhf.infer.intervals.upper_limits.upper_limit(
+        data, model, scan, level, return_results, **hypotest_kwargs
+    )
diff --git a/src/pyhf/infer/intervals.py → src/pyhf/infer/intervals/upper_limits.py b/src/pyhf/infer/intervals.py → src/pyhf/infer/intervals/upper_limits.py
@@ -7,7 +7,7 @@
 from pyhf import get_backend
 from pyhf.infer import hypotest
 
-__all__ = ["upperlimit", "upperlimit_auto", "upperlimit_fixed_scan"]
+__all__ = ["upper_limit", "linear_grid_scan", "toms748_scan"]
 
 
 def __dir__():
@@ -19,15 +19,15 @@ def _interp(x, xp, fp):
     return tb.astensor(np.interp(x, xp.tolist(), fp.tolist()))
 
 
-def upperlimit_auto(
+def toms748_scan(
     data,
     model,
-    low,
-    high,
+    bounds_low,
+    bounds_up,
     level=0.05,
     atol=2e-12,
     rtol=None,
-    from_upperlimit_fn=False,
+    from_upper_limit_fn=False,
     **hypotest_kwargs,
 ):
     """
@@ -44,7 +44,7 @@ def upperlimit_auto(
         ... )
         >>> observations = [51, 48]
         >>> data = pyhf.tensorlib.astensor(observations + model.config.auxdata)
-        >>> obs_limit, exp_limits = pyhf.infer.intervals.upperlimit_auto(
+        >>> obs_limit, exp_limits = pyhf.infer.intervals.upper_limits.toms748_scan(
         ...     data, model, 0., 5., rtol=0.01
         ... )
         >>> obs_limit
@@ -55,8 +55,8 @@ def upperlimit_auto(
     Args:
         data (:obj:`tensor`): The observed data.
         model (~pyhf.pdf.Model): The statistical model adhering to the schema ``model.json``.
-        low (:obj:`float`): Lower boundary of search region
-        high (:obj:`float`): Higher boundary of search region
+        bounds_low (:obj:`float`): Lower boundary of search interval.
+        bounds_up (:obj:`float`): Upper boundary of search interval.
         level (:obj:`float`): The threshold value to evaluate the interpolated results at.
                               Defaults to ``0.05``.
         atol (:obj:`float`): Absolute tolerance.
@@ -75,11 +75,14 @@ def upperlimit_auto(
 
             - Tensor: The observed upper limit on the POI.
             - Tensor: The expected upper limits on the POI.
+
+    .. versionadded:: 0.7.0
     """
     if rtol is None:
         rtol = 1e-15
         warn(
-            f"upperlimit_auto: rtol not provided, defaulting to {rtol}. For optimal performance rtol should be set to the highest acceptable relative tolerance."
+            f"toms748_scan: rtol not provided, defaulting to {rtol}.\n"
+            "For optimal performance rtol should be set to the highest acceptable relative tolerance."
         )
 
     cache = {}
@@ -120,35 +123,40 @@ def best_bracket(limit):
         upper = ks[neg][np.argmax(vals[neg])]
         return (lower, upper)
 
-    # extend low and high if they don't bracket CLs level
-    low_res = f_cached(low)
-    while np.any(np.array(low_res[0] + low_res[1]) < level):
-        low /= 2
-        low_res = f_cached(low)
-    high_res = f_cached(high)
-    while np.any(np.array(high_res[0] + high_res[1]) > level):
-        high *= 2
-        high_res = f_cached(high)
+    # extend bounds_low and bounds_up if they don't bracket CLs level
+    lower_results = f_cached(bounds_low)
+    # {lower,upper}_results[0] is an array and {lower,upper}_results[1] is a
+    # list of arrays so need to turn {lower,upper}_results[0] into list to
+    # concatenate them
+    while np.any(np.asarray([lower_results[0]] + lower_results[1]) < level):
+        bounds_low /= 2
+        lower_results = f_cached(bounds_low)
+    upper_results = f_cached(bounds_up)
+    while np.any(np.asarray([upper_results[0]] + upper_results[1]) > level):
+        bounds_up *= 2
+        upper_results = f_cached(bounds_up)
 
     tb, _ = get_backend()
-    obs = tb.astensor(toms748(f, low, high, args=(level, 0), k=2, xtol=atol, rtol=rtol))
+    obs = tb.astensor(
+        toms748(f, bounds_low, bounds_up, args=(level, 0), k=2, xtol=atol, rtol=rtol)
+    )
     exp = [
         tb.astensor(
-            toms748(f, *best_bracket(i), args=(level, i), k=2, xtol=atol, rtol=rtol)
+            toms748(f, *best_bracket(idx), args=(level, idx), k=2, xtol=atol, rtol=rtol)
         )
-        for i in range(1, 6)
+        for idx in range(1, 6)
     ]
-    if from_upperlimit_fn:
+    if from_upper_limit_fn:
         return obs, exp, (list(cache.keys()), list(cache.values()))
     return obs, exp
 
 
-def upperlimit_fixed_scan(
+def linear_grid_scan(
     data, model, scan, level=0.05, return_results=False, **hypotest_kwargs
 ):
     """
     Calculate an upper limit interval ``(0, poi_up)`` for a single
-    Parameter of Interest (POI) using a fixed scan through POI-space.
+    Parameter of Interest (POI) using a linear scan through POI-space.
 
     Example:
         >>> import numpy as np
@@ -160,7 +168,7 @@ def upperlimit_fixed_scan(
         >>> observations = [51, 48]
         >>> data = pyhf.tensorlib.astensor(observations + model.config.auxdata)
         >>> scan = np.linspace(0, 5, 21)
-        >>> obs_limit, exp_limits, (scan, results) = pyhf.infer.intervals.upperlimit(
+        >>> obs_limit, exp_limits, (scan, results) = pyhf.infer.intervals.upper_limits.upper_limit(
         ...     data, model, scan, return_results=True
         ... )
         >>> obs_limit
@@ -185,6 +193,8 @@ def upperlimit_fixed_scan(
             - Tuple of Tensors: The given ``scan`` along with the
               :class:`~pyhf.infer.hypotest` results at each test POI.
               Only returned when ``return_results`` is ``True``.
+
+    .. versionadded:: 0.7.0
     """
     tb, _ = get_backend()
     results = [
@@ -205,12 +215,12 @@ def upperlimit_fixed_scan(
     return obs_limit, exp_limits
 
 
-def upperlimit(
+def upper_limit(
     data, model, scan=None, level=0.05, return_results=False, **hypotest_kwargs
 ):
     """
     Calculate an upper limit interval ``(0, poi_up)`` for a single Parameter of
-    Interest (POI) using root-finding or a fixed scan through POI-space.
+    Interest (POI) using root-finding or a linear scan through POI-space.
 
     Example:
         >>> import numpy as np
@@ -222,7 +232,7 @@ def upperlimit(
         >>> observations = [51, 48]
         >>> data = pyhf.tensorlib.astensor(observations + model.config.auxdata)
         >>> scan = np.linspace(0, 5, 21)
-        >>> obs_limit, exp_limits, (scan, results) = pyhf.infer.intervals.upperlimit(
+        >>> obs_limit, exp_limits, (scan, results) = pyhf.infer.intervals.upper_limits.upper_limit(
         ...     data, model, scan, return_results=True
         ... )
         >>> obs_limit
@@ -233,7 +243,8 @@ def upperlimit(
     Args:
         data (:obj:`tensor`): The observed data.
         model (~pyhf.pdf.Model): The statistical model adhering to the schema ``model.json``.
-        scan (:obj:`iterable` or ``None``): Iterable of POI values or ``None`` to use ``upperlimit_auto``.
+        scan (:obj:`iterable` or ``None``): Iterable of POI values or ``None`` to use
+         :class:`~pyhf.infer.intervals.upper_limits.toms748_scan`.
         level (:obj:`float`): The threshold value to evaluate the interpolated results at.
         return_results (:obj:`bool`): Whether to return the per-point results.
 
@@ -245,22 +256,25 @@ def upperlimit(
             - Tuple of Tensors: The given ``scan`` along with the
               :class:`~pyhf.infer.hypotest` results at each test POI.
               Only returned when ``return_results`` is ``True``.
+
+    .. versionadded:: 0.7.0
     """
     if scan is not None:
-        return upperlimit_fixed_scan(
+        return linear_grid_scan(
             data, model, scan, level, return_results, **hypotest_kwargs
         )
     # else:
     bounds = model.config.suggested_bounds()[
         model.config.par_slice(model.config.poi_name).start
     ]
-    obs_limit, exp_limit, results = upperlimit_auto(
+    relative_tolerance = hypotest_kwargs.pop("rtol", 1e-8)
+    obs_limit, exp_limit, results = toms748_scan(
         data,
         model,
         bounds[0],
         bounds[1],
-        rtol=1e-3,
-        from_upperlimit_fn=True,
+        rtol=relative_tolerance,
+        from_upper_limit_fn=True,
         **hypotest_kwargs,
     )
     if return_results:

diff --git a/tests/constraints.txt b/tests/constraints.txt
@@ -1,5 +1,5 @@
 # core
-scipy==1.1.0
+scipy==1.2.0  # c.f. PR #1274
 click==8.0.0  # c.f. PR #1958, #1909
 tqdm==4.56.0
 jsonschema==4.15.0  # c.f. PR #1979

diff --git a/tests/contrib/baseline/test_plot_results_no_axis.png b/tests/contrib/baseline/test_plot_results_no_axis.png
diff --git a/tests/contrib/test_viz.py b/tests/contrib/test_viz.py
@@ -1,4 +1,5 @@
 import json
+import sys
 
 import matplotlib
 import matplotlib.pyplot as plt
@@ -67,6 +68,10 @@ def test_plot_results(datadir):
 
 
 @pytest.mark.mpl_image_compare
+@pytest.mark.xfail(
+    sys.version_info < (3, 8),
+    reason="baseline image generated with matplotlib v3.6.0 which is Python 3.8+",
+)
 def test_plot_results_no_axis(datadir):
     data = json.load(datadir.joinpath("hypotest_results.json").open(encoding="utf-8"))