We are currently solving the following issue within our repository. Here is the issue text:
--- BEGIN ISSUE ---
Bug in scipy minimization: init_pars and fixed_vals are not synced
# Description

Related: #1051 (found as part of this PR).

There is a bug in `pyhf`'s codebase (since ~forever) where the `init_par` for a fixed parameter might differ from the constrained value set for that fixed parameter. See for example v0.4.4:

https://github.com/scikit-hep/pyhf/blob/6ac0f6280c56111c32f46ec9aaad0286d05d36e7/src/pyhf/infer/mle.py#L52-L90
https://github.com/scikit-hep/pyhf/blob/6ac0f6280c56111c32f46ec9aaad0286d05d36e7/src/pyhf/optimize/opt_scipy.py#L16-L53

# Expected Behavior

I would've thought the underlying optimizers (particularly scipy) would be smart enough to handle or recognize when the constraint/init_pars are at odds with each other.

# Actual Behavior

See description.

# Steps to Reproduce

Follow along with me:

call `fixed_poi_fit` with:
- `poi_val = 0.5`
- `data`
- `pdf`
- `init_pars=[1, 1, 1]`

This will call `opt.minimize` with:

- `objective=twice_nll`
- `data=data`
- `pdf=pdf`
- `init_pars=[1, 1, 1]`
- `fixed_vals=[(0, 0.5)]`

at this point, `opt.minimize` (in the case of scipy) will do

```python
        constraints = [{'type': 'eq', 'fun': lambda v: v[indices] - values}]
```

which becomes essentially

```python
        constraints = [{'type': 'eq', 'fun': lambda v: v[0] - 0.5}]
```

but still calls `scipy.optimize.minimize` with `init_pars=[1, 1, 1]`. If one changes this to `[0.5, 1, 1]` -- one observes a slightly different result... but probably more correct.

Minuit does not have the same bug here because of https://github.com/scikit-hep/pyhf/blob/6ac0f6280c56111c32f46ec9aaad0286d05d36e7/src/pyhf/optimize/opt_minuit.py#L45 which updates the initval with the fixed value `initvals['p{}'.format(index)] = value`.

I was worried that the major refactor I did #951 broke this or changed this. In fact, I managed to keep the same bug in scipy, and did not introduce this bug into minuit.

--- END ISSUE ---

Below are some code segments, each from a relevant file. One or more of these files may contain bugs.
--- BEGIN FILES ---
Path: `src/pyhf/optimize/opt_scipy.py`
Content:
```
1 """SciPy Optimizer Class."""
2 from .. import exceptions
3 from .mixins import OptimizerMixin
4 import scipy
5 
6 
7 class scipy_optimizer(OptimizerMixin):
8     """
9     Optimizer that uses :func:`scipy.optimize.minimize`.
10     """
11 
12     __slots__ = ['name']
13 
14     def __init__(self, *args, **kwargs):
15         """
16         Initialize the scipy_optimizer.
17 
18         See :class:`pyhf.optimize.mixins.OptimizerMixin` for configuration options.
19         """
20         self.name = 'scipy'
21         super().__init__(*args, **kwargs)
22 
23     def _get_minimizer(
24         self, objective_and_grad, init_pars, init_bounds, fixed_vals=None, do_grad=False
25     ):
26         return scipy.optimize.minimize
27 
28     def _minimize(
29         self,
30         minimizer,
31         func,
32         x0,
33         do_grad=False,
34         bounds=None,
35         fixed_vals=None,
36         return_uncertainties=False,
37         options={},
38     ):
39         """
40         Same signature as :func:`scipy.optimize.minimize`.
41 
42         Minimizer Options:
43             maxiter (`int`): maximum number of iterations. Default is 100000.
44             verbose (`bool`): print verbose output during minimization. Default is off.
45             method (`str`): minimization routine. Default is 'SLSQP'.
46 
47         Returns:
48             fitresult (scipy.optimize.OptimizeResult): the fit result
49         """
50         maxiter = options.pop('maxiter', self.maxiter)
51         verbose = options.pop('verbose', self.verbose)
52         method = options.pop('method', 'SLSQP')
53         if options:
54             raise exceptions.Unsupported(
55                 f"Unsupported options were passed in: {list(options.keys())}."
56             )
57 
58         fixed_vals = fixed_vals or []
59         indices = [i for i, _ in fixed_vals]
60         values = [v for _, v in fixed_vals]
61         if fixed_vals:
62             constraints = [{'type': 'eq', 'fun': lambda v: v[indices] - values}]
63         else:
64             constraints = []
65 
66         return minimizer(
67             func,
68             x0,
69             method=method,
70             jac=do_grad,
71             bounds=bounds,
72             constraints=constraints,
73             options=dict(maxiter=maxiter, disp=bool(verbose)),
74         )
75
```

Path: `src/pyhf/infer/__init__.py`
Content:
```
1 """Inference for Statistical Models."""
2 
3 from .test_statistics import qmu
4 from .. import get_backend
5 from .calculators import AsymptoticCalculator
6 
7 
8 def hypotest(
9     poi_test, data, pdf, init_pars=None, par_bounds=None, qtilde=False, **kwargs
10 ):
11     r"""
12     Compute :math:`p`-values and test statistics for a single value of the parameter of interest.
13 
14     Example:
15         >>> import pyhf
16         >>> pyhf.set_backend("numpy")
17         >>> model = pyhf.simplemodels.hepdata_like(
18         ...     signal_data=[12.0, 11.0], bkg_data=[50.0, 52.0], bkg_uncerts=[3.0, 7.0]
19         ... )
20         >>> observations = [51, 48]
21         >>> data = pyhf.tensorlib.astensor(observations + model.config.auxdata)
22         >>> test_poi = 1.0
23         >>> CLs_obs, CLs_exp_band = pyhf.infer.hypotest(
24         ...     test_poi, data, model, qtilde=True, return_expected_set=True
25         ... )
26         >>> CLs_obs
27         array(0.05251554)
28         >>> CLs_exp_band
29         [array(0.00260641), array(0.01382066), array(0.06445521), array(0.23526104), array(0.57304182)]
30 
31     Args:
32         poi_test (Number or Tensor): The value of the parameter of interest (POI)
33         data (Number or Tensor): The data considered
34         pdf (~pyhf.pdf.Model): The statistical model adhering to the schema ``model.json``
35         init_pars (Array or Tensor): The initial parameter values to be used for minimization
36         par_bounds (Array or Tensor): The parameter value bounds to be used for minimization
37         qtilde (Bool): When ``True`` perform the calculation using the alternative
38          test statistic, :math:`\tilde{q}_{\mu}`, as defined under the Wald
39          approximation in Equation (62) of :xref:`arXiv:1007.1727`.
40 
41     Keyword Args:
42         return_tail_probs (bool): Bool for returning :math:`\mathrm{CL}_{s+b}` and :math:`\mathrm{CL}_{b}`
43         return_expected (bool): Bool for returning :math:`\mathrm{CL}_{\mathrm{exp}}`
44         return_expected_set (bool): Bool for returning the :math:`(-2,-1,0,1,2)\sigma` :math:`\mathrm{CL}_{\mathrm{exp}}` --- the "Brazil band"
45 
46     Returns:
47         Tuple of Floats and lists of Floats:
48 
49             - :math:`\mathrm{CL}_{s}`: The modified :math:`p`-value compared to
50               the given threshold :math:`\alpha`, typically taken to be :math:`0.05`,
51               defined in :xref:`arXiv:1007.1727` as
52 
53             .. math::
54 
55                 \mathrm{CL}_{s} = \frac{\mathrm{CL}_{s+b}}{\mathrm{CL}_{b}} = \frac{p_{s+b}}{1-p_{b}}
56 
57             to protect against excluding signal models in which there is little
58             sensitivity. In the case that :math:`\mathrm{CL}_{s} \leq \alpha`
59             the given signal model is excluded.
60 
61             - :math:`\left[\mathrm{CL}_{s+b}, \mathrm{CL}_{b}\right]`: The
62               signal + background model hypothesis :math:`p`-value
63 
64             .. math::
65 
66                 \mathrm{CL}_{s+b} = p_{s+b}
67                 = p\left(q \geq q_{\mathrm{obs}}\middle|s+b\right)
68                 = \int\limits_{q_{\mathrm{obs}}}^{\infty} f\left(q\,\middle|s+b\right)\,dq
69                 = 1 - F\left(q_{\mathrm{obs}}(\mu)\,\middle|\mu'\right)
70 
71             and 1 minus the background only model hypothesis :math:`p`-value
72 
73             .. math::
74 
75                 \mathrm{CL}_{b} = 1- p_{b}
76                 = p\left(q \geq q_{\mathrm{obs}}\middle|b\right)
77                 = 1 - \int\limits_{-\infty}^{q_{\mathrm{obs}}} f\left(q\,\middle|b\right)\,dq
78                 = 1 - F\left(q_{\mathrm{obs}}(\mu)\,\middle|0\right)
79 
80             for signal strength :math:`\mu` and model hypothesis signal strength
81             :math:`\mu'`, where the cumulative density functions
82             :math:`F\left(q(\mu)\,\middle|\mu'\right)` are given by Equations (57)
83             and (65) of :xref:`arXiv:1007.1727` for upper-limit-like test
84             statistic :math:`q \in \{q_{\mu}, \tilde{q}_{\mu}\}`.
85             Only returned when ``return_tail_probs`` is ``True``.
86 
87             .. note::
88 
89                 The definitions of the :math:`\mathrm{CL}_{s+b}` and
90                 :math:`\mathrm{CL}_{b}` used are based on profile likelihood
91                 ratio test statistics.
92                 This procedure is common in the LHC-era, but differs from
93                 procedures used in the LEP and Tevatron eras, as briefly
94                 discussed in :math:`\S` 3.8 of :xref:`arXiv:1007.1727`.
95 
96             - :math:`\mathrm{CL}_{s,\mathrm{exp}}`: The expected :math:`\mathrm{CL}_{s}`
97               value corresponding to the test statistic under the background
98               only hypothesis :math:`\left(\mu=0\right)`.
99               Only returned when ``return_expected`` is ``True``.
100 
101             - :math:`\mathrm{CL}_{s,\mathrm{exp}}` band: The set of expected
102               :math:`\mathrm{CL}_{s}` values corresponding to the median
103               significance of variations of the signal strength from the
104               background only hypothesis :math:`\left(\mu=0\right)` at
105               :math:`(-2,-1,0,1,2)\sigma`.
106               That is, the :math:`p`-values that satisfy Equation (89) of
107               :xref:`arXiv:1007.1727`
108 
109             .. math::
110 
111                 \mathrm{band}_{N\sigma} = \mu' + \sigma\,\Phi^{-1}\left(1-\alpha\right) \pm N\sigma
112 
113             for :math:`\mu'=0` and :math:`N \in \left\{-2, -1, 0, 1, 2\right\}`.
114             These values define the boundaries of an uncertainty band sometimes
115             referred to as the "Brazil band".
116             Only returned when ``return_expected_set`` is ``True``.
117 
118     """
119     init_pars = init_pars or pdf.config.suggested_init()
120     par_bounds = par_bounds or pdf.config.suggested_bounds()
121     tensorlib, _ = get_backend()
122 
123     calc = AsymptoticCalculator(data, pdf, init_pars, par_bounds, qtilde=qtilde)
124     teststat = calc.teststatistic(poi_test)
125     sig_plus_bkg_distribution, b_only_distribution = calc.distributions(poi_test)
126 
127     CLsb = sig_plus_bkg_distribution.pvalue(teststat)
128     CLb = b_only_distribution.pvalue(teststat)
129     CLs = CLsb / CLb
130     # Ensure that all CL values are 0-d tensors
131     CLsb, CLb, CLs = (
132         tensorlib.astensor(CLsb),
133         tensorlib.astensor(CLb),
134         tensorlib.astensor(CLs),
135     )
136 
137     _returns = [CLs]
138     if kwargs.get('return_tail_probs'):
139         _returns.append([CLsb, CLb])
140     if kwargs.get('return_expected_set'):
141         CLs_exp = []
142         for n_sigma in [2, 1, 0, -1, -2]:
143 
144             expected_bonly_teststat = b_only_distribution.expected_value(n_sigma)
145 
146             CLs = sig_plus_bkg_distribution.pvalue(
147                 expected_bonly_teststat
148             ) / b_only_distribution.pvalue(expected_bonly_teststat)
149             CLs_exp.append(tensorlib.astensor(CLs))
150         if kwargs.get('return_expected'):
151             _returns.append(CLs_exp[2])
152         _returns.append(CLs_exp)
153     elif kwargs.get('return_expected'):
154         n_sigma = 0
155         expected_bonly_teststat = b_only_distribution.expected_value(n_sigma)
156 
157         CLs = sig_plus_bkg_distribution.pvalue(
158             expected_bonly_teststat
159         ) / b_only_distribution.pvalue(expected_bonly_teststat)
160         _returns.append(tensorlib.astensor(CLs))
161     # Enforce a consistent return type of the observed CLs
162     return tuple(_returns) if len(_returns) > 1 else _returns[0]
163 
164 
165 __all__ = ['qmu', 'hypotest']
166
```
--- END FILES ---

Please first localize the bug based on the issue statement, and then generate a patch according to the `git diff` format fenced by three backticks.

Here is an example:

```diff
diff --git a/examples/server_async.py b/examples/server_async.py
--- a/examples/server_async.py
+++ b/examples/server_async.py
@@ -313,4 +313,4 @@
 
 
 if __name__ == "__main__":
-    asyncio.run(run_async_server("."), debug=True)
+    asyncio.run(run_async_server(), debug=True)
diff --git a/examples/server_sync.py b/examples/server_sync.py
--- a/examples/server_sync.py
+++ b/examples/server_sync.py
@@ -313,5 +313,5 @@
 
 
 if __name__ == "__main__":
-    server = run_sync_server(".")
+    server = run_sync_server()
     server.shutdown()

```

diff --git a/src/pyhf/infer/__init__.py b/src/pyhf/infer/__init__.py
--- a/src/pyhf/infer/__init__.py
+++ b/src/pyhf/infer/__init__.py
@@ -24,9 +24,9 @@
         ...     test_poi, data, model, qtilde=True, return_expected_set=True
         ... )
         >>> CLs_obs
-        array(0.05251554)
+        array(0.05251497)
         >>> CLs_exp_band
-        [array(0.00260641), array(0.01382066), array(0.06445521), array(0.23526104), array(0.57304182)]
+        [array(0.00260626), array(0.01382005), array(0.06445321), array(0.23525644), array(0.57303621)]
 
     Args:
         poi_test (Number or Tensor): The value of the parameter of interest (POI)
diff --git a/src/pyhf/optimize/opt_scipy.py b/src/pyhf/optimize/opt_scipy.py
--- a/src/pyhf/optimize/opt_scipy.py
+++ b/src/pyhf/optimize/opt_scipy.py
@@ -60,6 +60,9 @@
         values = [v for _, v in fixed_vals]
         if fixed_vals:
             constraints = [{'type': 'eq', 'fun': lambda v: v[indices] - values}]
+            # update the initial values to the fixed value for any fixed parameter
+            for idx, fixed_val in fixed_vals:
+                x0[idx] = fixed_val
         else:
             constraints = []

{"golden_diff": "diff --git a/src/pyhf/infer/__init__.py b/src/pyhf/infer/__init__.py\n--- a/src/pyhf/infer/__init__.py\n+++ b/src/pyhf/infer/__init__.py\n@@ -24,9 +24,9 @@\n         ...     test_poi, data, model, qtilde=True, return_expected_set=True\n         ... )\n         >>> CLs_obs\n-        array(0.05251554)\n+        array(0.05251497)\n         >>> CLs_exp_band\n-        [array(0.00260641), array(0.01382066), array(0.06445521), array(0.23526104), array(0.57304182)]\n+        [array(0.00260626), array(0.01382005), array(0.06445321), array(0.23525644), array(0.57303621)]\n \n     Args:\n         poi_test (Number or Tensor): The value of the parameter of interest (POI)\ndiff --git a/src/pyhf/optimize/opt_scipy.py b/src/pyhf/optimize/opt_scipy.py\n--- a/src/pyhf/optimize/opt_scipy.py\n+++ b/src/pyhf/optimize/opt_scipy.py\n@@ -60,6 +60,9 @@\n         values = [v for _, v in fixed_vals]\n         if fixed_vals:\n             constraints = [{'type': 'eq', 'fun': lambda v: v[indices] - values}]\n+            # update the initial values to the fixed value for any fixed parameter\n+            for idx, fixed_val in fixed_vals:\n+                x0[idx] = fixed_val\n         else:\n             constraints = []\n", "issue": "Bug in scipy minimization: init_pars and fixed_vals are not synced\n# Description\r\n\r\nRelated: #1051 (found as part of this PR).\r\n\r\nThere is a bug in `pyhf`'s codebase (since ~forever) where the `init_par` for a fixed parameter might differ from the constrained value set for that fixed parameter. See for example v0.4.4:\r\n\r\nhttps://github.com/scikit-hep/pyhf/blob/6ac0f6280c56111c32f46ec9aaad0286d05d36e7/src/pyhf/infer/mle.py#L52-L90\r\nhttps://github.com/scikit-hep/pyhf/blob/6ac0f6280c56111c32f46ec9aaad0286d05d36e7/src/pyhf/optimize/opt_scipy.py#L16-L53\r\n\r\n# Expected Behavior\r\n\r\nI would've thought the underlying optimizers (particularly scipy) would be smart enough to handle or recognize when the constraint/init_pars are at odds with each other.\r\n\r\n# Actual Behavior\r\n\r\nSee description.\r\n\r\n# Steps to Reproduce\r\n\r\nFollow along with me:\r\n\r\ncall `fixed_poi_fit` with:\r\n- `poi_val = 0.5`\r\n- `data`\r\n- `pdf`\r\n- `init_pars=[1, 1, 1]`\r\n\r\nThis will call `opt.minimize` with:\r\n\r\n- `objective=twice_nll`\r\n- `data=data`\r\n- `pdf=pdf`\r\n- `init_pars=[1, 1, 1]`\r\n- `fixed_vals=[(0, 0.5)]`\r\n\r\nat this point, `opt.minimize` (in the case of scipy) will do\r\n\r\n```python\r\n        constraints = [{'type': 'eq', 'fun': lambda v: v[indices] - values}]\r\n```\r\n\r\nwhich becomes essentially\r\n\r\n```python\r\n        constraints = [{'type': 'eq', 'fun': lambda v: v[0] - 0.5}]\r\n```\r\n\r\nbut still calls `scipy.optimize.minimize` with `init_pars=[1, 1, 1]`. If one changes this to `[0.5, 1, 1]` -- one observes a slightly different result... but probably more correct.\r\n\r\nMinuit does not have the same bug here because of https://github.com/scikit-hep/pyhf/blob/6ac0f6280c56111c32f46ec9aaad0286d05d36e7/src/pyhf/optimize/opt_minuit.py#L45 which updates the initval with the fixed value `initvals['p{}'.format(index)] = value`.\r\n\r\nI was worried that the major refactor I did #951 broke this or changed this. In fact, I managed to keep the same bug in scipy, and did not introduce this bug into minuit.\n", "before_files": [{"content": "\"\"\"SciPy Optimizer Class.\"\"\"\nfrom .. import exceptions\nfrom .mixins import OptimizerMixin\nimport scipy\n\n\nclass scipy_optimizer(OptimizerMixin):\n    \"\"\"\n    Optimizer that uses :func:`scipy.optimize.minimize`.\n    \"\"\"\n\n    __slots__ = ['name']\n\n    def __init__(self, *args, **kwargs):\n        \"\"\"\n        Initialize the scipy_optimizer.\n\n        See :class:`pyhf.optimize.mixins.OptimizerMixin` for configuration options.\n        \"\"\"\n        self.name = 'scipy'\n        super().__init__(*args, **kwargs)\n\n    def _get_minimizer(\n        self, objective_and_grad, init_pars, init_bounds, fixed_vals=None, do_grad=False\n    ):\n        return scipy.optimize.minimize\n\n    def _minimize(\n        self,\n        minimizer,\n        func,\n        x0,\n        do_grad=False,\n        bounds=None,\n        fixed_vals=None,\n        return_uncertainties=False,\n        options={},\n    ):\n        \"\"\"\n        Same signature as :func:`scipy.optimize.minimize`.\n\n        Minimizer Options:\n            maxiter (`int`): maximum number of iterations. Default is 100000.\n            verbose (`bool`): print verbose output during minimization. Default is off.\n            method (`str`): minimization routine. Default is 'SLSQP'.\n\n        Returns:\n            fitresult (scipy.optimize.OptimizeResult): the fit result\n        \"\"\"\n        maxiter = options.pop('maxiter', self.maxiter)\n        verbose = options.pop('verbose', self.verbose)\n        method = options.pop('method', 'SLSQP')\n        if options:\n            raise exceptions.Unsupported(\n                f\"Unsupported options were passed in: {list(options.keys())}.\"\n            )\n\n        fixed_vals = fixed_vals or []\n        indices = [i for i, _ in fixed_vals]\n        values = [v for _, v in fixed_vals]\n        if fixed_vals:\n            constraints = [{'type': 'eq', 'fun': lambda v: v[indices] - values}]\n        else:\n            constraints = []\n\n        return minimizer(\n            func,\n            x0,\n            method=method,\n            jac=do_grad,\n            bounds=bounds,\n            constraints=constraints,\n            options=dict(maxiter=maxiter, disp=bool(verbose)),\n        )\n", "path": "src/pyhf/optimize/opt_scipy.py"}, {"content": "\"\"\"Inference for Statistical Models.\"\"\"\n\nfrom .test_statistics import qmu\nfrom .. import get_backend\nfrom .calculators import AsymptoticCalculator\n\n\ndef hypotest(\n    poi_test, data, pdf, init_pars=None, par_bounds=None, qtilde=False, **kwargs\n):\n    r\"\"\"\n    Compute :math:`p`-values and test statistics for a single value of the parameter of interest.\n\n    Example:\n        >>> import pyhf\n        >>> pyhf.set_backend(\"numpy\")\n        >>> model = pyhf.simplemodels.hepdata_like(\n        ...     signal_data=[12.0, 11.0], bkg_data=[50.0, 52.0], bkg_uncerts=[3.0, 7.0]\n        ... )\n        >>> observations = [51, 48]\n        >>> data = pyhf.tensorlib.astensor(observations + model.config.auxdata)\n        >>> test_poi = 1.0\n        >>> CLs_obs, CLs_exp_band = pyhf.infer.hypotest(\n        ...     test_poi, data, model, qtilde=True, return_expected_set=True\n        ... )\n        >>> CLs_obs\n        array(0.05251554)\n        >>> CLs_exp_band\n        [array(0.00260641), array(0.01382066), array(0.06445521), array(0.23526104), array(0.57304182)]\n\n    Args:\n        poi_test (Number or Tensor): The value of the parameter of interest (POI)\n        data (Number or Tensor): The data considered\n        pdf (~pyhf.pdf.Model): The statistical model adhering to the schema ``model.json``\n        init_pars (Array or Tensor): The initial parameter values to be used for minimization\n        par_bounds (Array or Tensor): The parameter value bounds to be used for minimization\n        qtilde (Bool): When ``True`` perform the calculation using the alternative\n         test statistic, :math:`\\tilde{q}_{\\mu}`, as defined under the Wald\n         approximation in Equation (62) of :xref:`arXiv:1007.1727`.\n\n    Keyword Args:\n        return_tail_probs (bool): Bool for returning :math:`\\mathrm{CL}_{s+b}` and :math:`\\mathrm{CL}_{b}`\n        return_expected (bool): Bool for returning :math:`\\mathrm{CL}_{\\mathrm{exp}}`\n        return_expected_set (bool): Bool for returning the :math:`(-2,-1,0,1,2)\\sigma` :math:`\\mathrm{CL}_{\\mathrm{exp}}` --- the \"Brazil band\"\n\n    Returns:\n        Tuple of Floats and lists of Floats:\n\n            - :math:`\\mathrm{CL}_{s}`: The modified :math:`p`-value compared to\n              the given threshold :math:`\\alpha`, typically taken to be :math:`0.05`,\n              defined in :xref:`arXiv:1007.1727` as\n\n            .. math::\n\n                \\mathrm{CL}_{s} = \\frac{\\mathrm{CL}_{s+b}}{\\mathrm{CL}_{b}} = \\frac{p_{s+b}}{1-p_{b}}\n\n            to protect against excluding signal models in which there is little\n            sensitivity. In the case that :math:`\\mathrm{CL}_{s} \\leq \\alpha`\n            the given signal model is excluded.\n\n            - :math:`\\left[\\mathrm{CL}_{s+b}, \\mathrm{CL}_{b}\\right]`: The\n              signal + background model hypothesis :math:`p`-value\n\n            .. math::\n\n                \\mathrm{CL}_{s+b} = p_{s+b}\n                = p\\left(q \\geq q_{\\mathrm{obs}}\\middle|s+b\\right)\n                = \\int\\limits_{q_{\\mathrm{obs}}}^{\\infty} f\\left(q\\,\\middle|s+b\\right)\\,dq\n                = 1 - F\\left(q_{\\mathrm{obs}}(\\mu)\\,\\middle|\\mu'\\right)\n\n            and 1 minus the background only model hypothesis :math:`p`-value\n\n            .. math::\n\n                \\mathrm{CL}_{b} = 1- p_{b}\n                = p\\left(q \\geq q_{\\mathrm{obs}}\\middle|b\\right)\n                = 1 - \\int\\limits_{-\\infty}^{q_{\\mathrm{obs}}} f\\left(q\\,\\middle|b\\right)\\,dq\n                = 1 - F\\left(q_{\\mathrm{obs}}(\\mu)\\,\\middle|0\\right)\n\n            for signal strength :math:`\\mu` and model hypothesis signal strength\n            :math:`\\mu'`, where the cumulative density functions\n            :math:`F\\left(q(\\mu)\\,\\middle|\\mu'\\right)` are given by Equations (57)\n            and (65) of :xref:`arXiv:1007.1727` for upper-limit-like test\n            statistic :math:`q \\in \\{q_{\\mu}, \\tilde{q}_{\\mu}\\}`.\n            Only returned when ``return_tail_probs`` is ``True``.\n\n            .. note::\n\n                The definitions of the :math:`\\mathrm{CL}_{s+b}` and\n                :math:`\\mathrm{CL}_{b}` used are based on profile likelihood\n                ratio test statistics.\n                This procedure is common in the LHC-era, but differs from\n                procedures used in the LEP and Tevatron eras, as briefly\n                discussed in :math:`\\S` 3.8 of :xref:`arXiv:1007.1727`.\n\n            - :math:`\\mathrm{CL}_{s,\\mathrm{exp}}`: The expected :math:`\\mathrm{CL}_{s}`\n              value corresponding to the test statistic under the background\n              only hypothesis :math:`\\left(\\mu=0\\right)`.\n              Only returned when ``return_expected`` is ``True``.\n\n            - :math:`\\mathrm{CL}_{s,\\mathrm{exp}}` band: The set of expected\n              :math:`\\mathrm{CL}_{s}` values corresponding to the median\n              significance of variations of the signal strength from the\n              background only hypothesis :math:`\\left(\\mu=0\\right)` at\n              :math:`(-2,-1,0,1,2)\\sigma`.\n              That is, the :math:`p`-values that satisfy Equation (89) of\n              :xref:`arXiv:1007.1727`\n\n            .. math::\n\n                \\mathrm{band}_{N\\sigma} = \\mu' + \\sigma\\,\\Phi^{-1}\\left(1-\\alpha\\right) \\pm N\\sigma\n\n            for :math:`\\mu'=0` and :math:`N \\in \\left\\{-2, -1, 0, 1, 2\\right\\}`.\n            These values define the boundaries of an uncertainty band sometimes\n            referred to as the \"Brazil band\".\n            Only returned when ``return_expected_set`` is ``True``.\n\n    \"\"\"\n    init_pars = init_pars or pdf.config.suggested_init()\n    par_bounds = par_bounds or pdf.config.suggested_bounds()\n    tensorlib, _ = get_backend()\n\n    calc = AsymptoticCalculator(data, pdf, init_pars, par_bounds, qtilde=qtilde)\n    teststat = calc.teststatistic(poi_test)\n    sig_plus_bkg_distribution, b_only_distribution = calc.distributions(poi_test)\n\n    CLsb = sig_plus_bkg_distribution.pvalue(teststat)\n    CLb = b_only_distribution.pvalue(teststat)\n    CLs = CLsb / CLb\n    # Ensure that all CL values are 0-d tensors\n    CLsb, CLb, CLs = (\n        tensorlib.astensor(CLsb),\n        tensorlib.astensor(CLb),\n        tensorlib.astensor(CLs),\n    )\n\n    _returns = [CLs]\n    if kwargs.get('return_tail_probs'):\n        _returns.append([CLsb, CLb])\n    if kwargs.get('return_expected_set'):\n        CLs_exp = []\n        for n_sigma in [2, 1, 0, -1, -2]:\n\n            expected_bonly_teststat = b_only_distribution.expected_value(n_sigma)\n\n            CLs = sig_plus_bkg_distribution.pvalue(\n                expected_bonly_teststat\n            ) / b_only_distribution.pvalue(expected_bonly_teststat)\n            CLs_exp.append(tensorlib.astensor(CLs))\n        if kwargs.get('return_expected'):\n            _returns.append(CLs_exp[2])\n        _returns.append(CLs_exp)\n    elif kwargs.get('return_expected'):\n        n_sigma = 0\n        expected_bonly_teststat = b_only_distribution.expected_value(n_sigma)\n\n        CLs = sig_plus_bkg_distribution.pvalue(\n            expected_bonly_teststat\n        ) / b_only_distribution.pvalue(expected_bonly_teststat)\n        _returns.append(tensorlib.astensor(CLs))\n    # Enforce a consistent return type of the observed CLs\n    return tuple(_returns) if len(_returns) > 1 else _returns[0]\n\n\n__all__ = ['qmu', 'hypotest']\n", "path": "src/pyhf/infer/__init__.py"}], "after_files": [{"content": "\"\"\"SciPy Optimizer Class.\"\"\"\nfrom .. import exceptions\nfrom .mixins import OptimizerMixin\nimport scipy\n\n\nclass scipy_optimizer(OptimizerMixin):\n    \"\"\"\n    Optimizer that uses :func:`scipy.optimize.minimize`.\n    \"\"\"\n\n    __slots__ = ['name']\n\n    def __init__(self, *args, **kwargs):\n        \"\"\"\n        Initialize the scipy_optimizer.\n\n        See :class:`pyhf.optimize.mixins.OptimizerMixin` for configuration options.\n        \"\"\"\n        self.name = 'scipy'\n        super().__init__(*args, **kwargs)\n\n    def _get_minimizer(\n        self, objective_and_grad, init_pars, init_bounds, fixed_vals=None, do_grad=False\n    ):\n        return scipy.optimize.minimize\n\n    def _minimize(\n        self,\n        minimizer,\n        func,\n        x0,\n        do_grad=False,\n        bounds=None,\n        fixed_vals=None,\n        return_uncertainties=False,\n        options={},\n    ):\n        \"\"\"\n        Same signature as :func:`scipy.optimize.minimize`.\n\n        Minimizer Options:\n            maxiter (`int`): maximum number of iterations. Default is 100000.\n            verbose (`bool`): print verbose output during minimization. Default is off.\n            method (`str`): minimization routine. Default is 'SLSQP'.\n\n        Returns:\n            fitresult (scipy.optimize.OptimizeResult): the fit result\n        \"\"\"\n        maxiter = options.pop('maxiter', self.maxiter)\n        verbose = options.pop('verbose', self.verbose)\n        method = options.pop('method', 'SLSQP')\n        if options:\n            raise exceptions.Unsupported(\n                f\"Unsupported options were passed in: {list(options.keys())}.\"\n            )\n\n        fixed_vals = fixed_vals or []\n        indices = [i for i, _ in fixed_vals]\n        values = [v for _, v in fixed_vals]\n        if fixed_vals:\n            constraints = [{'type': 'eq', 'fun': lambda v: v[indices] - values}]\n            # update the initial values to the fixed value for any fixed parameter\n            for idx, fixed_val in fixed_vals:\n                x0[idx] = fixed_val\n        else:\n            constraints = []\n\n        return minimizer(\n            func,\n            x0,\n            method=method,\n            jac=do_grad,\n            bounds=bounds,\n            constraints=constraints,\n            options=dict(maxiter=maxiter, disp=bool(verbose)),\n        )\n", "path": "src/pyhf/optimize/opt_scipy.py"}, {"content": "\"\"\"Inference for Statistical Models.\"\"\"\n\nfrom .test_statistics import qmu\nfrom .. import get_backend\nfrom .calculators import AsymptoticCalculator\n\n\ndef hypotest(\n    poi_test, data, pdf, init_pars=None, par_bounds=None, qtilde=False, **kwargs\n):\n    r\"\"\"\n    Compute :math:`p`-values and test statistics for a single value of the parameter of interest.\n\n    Example:\n        >>> import pyhf\n        >>> pyhf.set_backend(\"numpy\")\n        >>> model = pyhf.simplemodels.hepdata_like(\n        ...     signal_data=[12.0, 11.0], bkg_data=[50.0, 52.0], bkg_uncerts=[3.0, 7.0]\n        ... )\n        >>> observations = [51, 48]\n        >>> data = pyhf.tensorlib.astensor(observations + model.config.auxdata)\n        >>> test_poi = 1.0\n        >>> CLs_obs, CLs_exp_band = pyhf.infer.hypotest(\n        ...     test_poi, data, model, qtilde=True, return_expected_set=True\n        ... )\n        >>> CLs_obs\n        array(0.05251497)\n        >>> CLs_exp_band\n        [array(0.00260626), array(0.01382005), array(0.06445321), array(0.23525644), array(0.57303621)]\n\n    Args:\n        poi_test (Number or Tensor): The value of the parameter of interest (POI)\n        data (Number or Tensor): The data considered\n        pdf (~pyhf.pdf.Model): The statistical model adhering to the schema ``model.json``\n        init_pars (Array or Tensor): The initial parameter values to be used for minimization\n        par_bounds (Array or Tensor): The parameter value bounds to be used for minimization\n        qtilde (Bool): When ``True`` perform the calculation using the alternative\n         test statistic, :math:`\\tilde{q}_{\\mu}`, as defined under the Wald\n         approximation in Equation (62) of :xref:`arXiv:1007.1727`.\n\n    Keyword Args:\n        return_tail_probs (bool): Bool for returning :math:`\\mathrm{CL}_{s+b}` and :math:`\\mathrm{CL}_{b}`\n        return_expected (bool): Bool for returning :math:`\\mathrm{CL}_{\\mathrm{exp}}`\n        return_expected_set (bool): Bool for returning the :math:`(-2,-1,0,1,2)\\sigma` :math:`\\mathrm{CL}_{\\mathrm{exp}}` --- the \"Brazil band\"\n\n    Returns:\n        Tuple of Floats and lists of Floats:\n\n            - :math:`\\mathrm{CL}_{s}`: The modified :math:`p`-value compared to\n              the given threshold :math:`\\alpha`, typically taken to be :math:`0.05`,\n              defined in :xref:`arXiv:1007.1727` as\n\n            .. math::\n\n                \\mathrm{CL}_{s} = \\frac{\\mathrm{CL}_{s+b}}{\\mathrm{CL}_{b}} = \\frac{p_{s+b}}{1-p_{b}}\n\n            to protect against excluding signal models in which there is little\n            sensitivity. In the case that :math:`\\mathrm{CL}_{s} \\leq \\alpha`\n            the given signal model is excluded.\n\n            - :math:`\\left[\\mathrm{CL}_{s+b}, \\mathrm{CL}_{b}\\right]`: The\n              signal + background model hypothesis :math:`p`-value\n\n            .. math::\n\n                \\mathrm{CL}_{s+b} = p_{s+b}\n                = p\\left(q \\geq q_{\\mathrm{obs}}\\middle|s+b\\right)\n                = \\int\\limits_{q_{\\mathrm{obs}}}^{\\infty} f\\left(q\\,\\middle|s+b\\right)\\,dq\n                = 1 - F\\left(q_{\\mathrm{obs}}(\\mu)\\,\\middle|\\mu'\\right)\n\n            and 1 minus the background only model hypothesis :math:`p`-value\n\n            .. math::\n\n                \\mathrm{CL}_{b} = 1- p_{b}\n                = p\\left(q \\geq q_{\\mathrm{obs}}\\middle|b\\right)\n                = 1 - \\int\\limits_{-\\infty}^{q_{\\mathrm{obs}}} f\\left(q\\,\\middle|b\\right)\\,dq\n                = 1 - F\\left(q_{\\mathrm{obs}}(\\mu)\\,\\middle|0\\right)\n\n            for signal strength :math:`\\mu` and model hypothesis signal strength\n            :math:`\\mu'`, where the cumulative density functions\n            :math:`F\\left(q(\\mu)\\,\\middle|\\mu'\\right)` are given by Equations (57)\n            and (65) of :xref:`arXiv:1007.1727` for upper-limit-like test\n            statistic :math:`q \\in \\{q_{\\mu}, \\tilde{q}_{\\mu}\\}`.\n            Only returned when ``return_tail_probs`` is ``True``.\n\n            .. note::\n\n                The definitions of the :math:`\\mathrm{CL}_{s+b}` and\n                :math:`\\mathrm{CL}_{b}` used are based on profile likelihood\n                ratio test statistics.\n                This procedure is common in the LHC-era, but differs from\n                procedures used in the LEP and Tevatron eras, as briefly\n                discussed in :math:`\\S` 3.8 of :xref:`arXiv:1007.1727`.\n\n            - :math:`\\mathrm{CL}_{s,\\mathrm{exp}}`: The expected :math:`\\mathrm{CL}_{s}`\n              value corresponding to the test statistic under the background\n              only hypothesis :math:`\\left(\\mu=0\\right)`.\n              Only returned when ``return_expected`` is ``True``.\n\n            - :math:`\\mathrm{CL}_{s,\\mathrm{exp}}` band: The set of expected\n              :math:`\\mathrm{CL}_{s}` values corresponding to the median\n              significance of variations of the signal strength from the\n              background only hypothesis :math:`\\left(\\mu=0\\right)` at\n              :math:`(-2,-1,0,1,2)\\sigma`.\n              That is, the :math:`p`-values that satisfy Equation (89) of\n              :xref:`arXiv:1007.1727`\n\n            .. math::\n\n                \\mathrm{band}_{N\\sigma} = \\mu' + \\sigma\\,\\Phi^{-1}\\left(1-\\alpha\\right) \\pm N\\sigma\n\n            for :math:`\\mu'=0` and :math:`N \\in \\left\\{-2, -1, 0, 1, 2\\right\\}`.\n            These values define the boundaries of an uncertainty band sometimes\n            referred to as the \"Brazil band\".\n            Only returned when ``return_expected_set`` is ``True``.\n\n    \"\"\"\n    init_pars = init_pars or pdf.config.suggested_init()\n    par_bounds = par_bounds or pdf.config.suggested_bounds()\n    tensorlib, _ = get_backend()\n\n    calc = AsymptoticCalculator(data, pdf, init_pars, par_bounds, qtilde=qtilde)\n    teststat = calc.teststatistic(poi_test)\n    sig_plus_bkg_distribution, b_only_distribution = calc.distributions(poi_test)\n\n    CLsb = sig_plus_bkg_distribution.pvalue(teststat)\n    CLb = b_only_distribution.pvalue(teststat)\n    CLs = CLsb / CLb\n    # Ensure that all CL values are 0-d tensors\n    CLsb, CLb, CLs = (\n        tensorlib.astensor(CLsb),\n        tensorlib.astensor(CLb),\n        tensorlib.astensor(CLs),\n    )\n\n    _returns = [CLs]\n    if kwargs.get('return_tail_probs'):\n        _returns.append([CLsb, CLb])\n    if kwargs.get('return_expected_set'):\n        CLs_exp = []\n        for n_sigma in [2, 1, 0, -1, -2]:\n\n            expected_bonly_teststat = b_only_distribution.expected_value(n_sigma)\n\n            CLs = sig_plus_bkg_distribution.pvalue(\n                expected_bonly_teststat\n            ) / b_only_distribution.pvalue(expected_bonly_teststat)\n            CLs_exp.append(tensorlib.astensor(CLs))\n        if kwargs.get('return_expected'):\n            _returns.append(CLs_exp[2])\n        _returns.append(CLs_exp)\n    elif kwargs.get('return_expected'):\n        n_sigma = 0\n        expected_bonly_teststat = b_only_distribution.expected_value(n_sigma)\n\n        CLs = sig_plus_bkg_distribution.pvalue(\n            expected_bonly_teststat\n        ) / b_only_distribution.pvalue(expected_bonly_teststat)\n        _returns.append(tensorlib.astensor(CLs))\n    # Enforce a consistent return type of the observed CLs\n    return tuple(_returns) if len(_returns) > 1 else _returns[0]\n\n\n__all__ = ['qmu', 'hypotest']\n", "path": "src/pyhf/infer/__init__.py"}]}