fix mixed-monotonic bug

Author: Joe-Vincent

binomial_cis/volume.py

@@ -163,9 +163,9 @@ def max_expected_shortage(alpha, n, tol=1e-3, verbose=True, randomized=True):
     """
     I = Interval(0,1)
     if randomized:
-        def F(p1, p2): expected_shortage_mixed_monotonic(llc_accept_prob, alpha, n, p1, p2)
+        def F(p1, p2): return expected_shortage_mixed_monotonic(llc_accept_prob, alpha, n, p1, p2)
     else:
-        def F(p1, p2): expected_shortage_mixed_monotonic_cp(alpha, n, p1, p2)
+        def F(p1, p2): return expected_shortage_mixed_monotonic_cp(alpha, n, p1, p2)
 
     ub, lb, p_lb, num_iters = mmp_solve(F, I, tol=tol, max_iters=1000, verbose=verbose)
     return ub, lb, p_lb, num_iters
@@ -224,7 +224,7 @@ def max_expected_width(alpha, n, tol=1e-3, verbose=True):
     I = Interval(0,1)
     # expected width is increasing in p2 and decreasing in p1
     # mmp_solve expects increasing in first arg and decreasing in second arg
-    def F(p2, p1): expected_width_mixed_monotonic(llc_accept_prob_2_sided, alpha, n, p1, p2)
+    def F(p2, p1): return expected_width_mixed_monotonic(llc_accept_prob_2_sided, alpha, n, p1, p2)
     ub, lb, p_lb, num_iters = mmp_solve(F, I, tol=tol, max_iters=1000, verbose=verbose)
     return ub, lb, p_lb, num_iters
 

docs/_include/tests.rst

@@ -20,9 +20,10 @@ Then, create a virtual environment and load the dependencies (these commands are
 
 binom_ci_test.py
 ******************
-This file contains functions to test the correctness of the lower, upper, and 2-sided bounds.
+This file contains functions to test the correctness of the lower, upper, 2-sided bounds, and the mixed monotonic solver.
 The lower/upper bound tests ensure the bounds we compute are better than Clopper-Pearson bounds (without being too optimistic).
 The 2-sided test ensures that our 2-sided bounds agree with those computed by Blyth and Hutchinson in their 1960 paper titled "Table of Neyman-Shortest Unbiased Confidence Intervals for the Binomial Parameter".
+The mixed monotonic test ensures that the computed maximum of expected shortage/excess/width is greater than the sample-based maximum of these quantities.
 We use a Github Action to automatically run the tests in this file via ``pytest``. 
 To run the tests yourself simply navigate to the ``tests/`` directory and run
 

docs/conf.py

@@ -13,7 +13,7 @@
 project = 'binomial_cis'
 copyright = '2024, Joe Vincent'
 author = 'Joe Vincent'
-release = '0.0.12'
+release = '0.0.11'
 
 # -- General configuration ---------------------------------------------------
 # https://www.sphinx-doc.org/en/master/usage/configuration.html#general-configuration

setup.py

@@ -2,7 +2,7 @@
 
 setup(
    name="binomial_cis",
-   version='0.0.12',
+   version='0.0.11',
    author="Joe Vincent",
    description="Confidence intervals for binomial distributions.",
    packages=find_packages(),

tests/README.md

@@ -1,6 +1,6 @@
 The tests are written in the form of interactive notebooks.
 
 - `2_side_validation.ipynb`: This notebook can be used to more easily inspect any differences between our 2-sided bounds and those from the Blyth paper.
-- `binom_ci_test.py`: This file can be used to test the correctness of the lower, upper, and 2-sided bounds. 
+- `binom_ci_test.py`: This file can be used to test the correctness of the lower, upper, and 2-sided bounds, and the correctness of the mixed-monotonic solver.
 - `binom_helper_validation.ipynb`: This notebook can be used to validate the binomial helper functions from `binomial_cis/binomial_helper.py`. Accuracy and speed of our implementation is tested against the SciPy implementation.
 - `conf_set_validation.ipynb`: This notebook can be used to validate the probabilistic guarantees of the confidence intervals using Monte Carlo simulation.

tests/binom_ci_test.py

@@ -5,6 +5,11 @@
 from math import isclose
 import pandas as pd
 
+from binomial_cis.conf_intervals import llc_accept_prob, llc_accept_prob_2_sided
+from binomial_cis.volume import expected_excess, max_expected_excess
+from binomial_cis.volume import expected_shortage, max_expected_shortage
+from binomial_cis.volume import expected_width, max_expected_width
+
 
 # define range of test conditions
 ns = np.array(range(1,50+1))
@@ -208,5 +213,22 @@ def test_2_sided():
 
 
 
-
+def test_mixed_monotonicity():
+    """
+    Test that our values for MES, MEE, and MEW actually upper bound the
+    sample-based maximum of these quantities
+    """
+    n = 10
+    alpha = 0.05
+    mes_ub, _, _, _ = max_expected_shortage(alpha, n, verbose=False)
+    mee_ub, _, _, _ = max_expected_excess(alpha, n, verbose=False)
+    mew_ub, _, _, _ = max_expected_width(alpha, n, verbose=False)
+
+    sample_mes = max([expected_shortage(llc_accept_prob, alpha, n, p) for p in np.linspace(0.01, 0.99, num=50)])
+    sample_mee = max([expected_excess(llc_accept_prob, alpha, n, p) for p in np.linspace(0.01, 0.99, num=50)])
+    sample_mew = max([expected_width(llc_accept_prob_2_sided, alpha, n, p) for p in np.linspace(0.01, 0.99, num=50)])
+
+    assert mes_ub >= sample_mes
+    assert mee_ub >= sample_mee
+    assert mew_ub >= sample_mew