From 35cfd28357dc70b1a2e32b1faf196d5c4e1a3348 Mon Sep 17 00:00:00 2001 From: Hans Dembinski Date: Mon, 13 Feb 2023 23:40:29 +0100 Subject: [PATCH] Fix coverage, typing, Template for 2D models (#836) --- src/iminuit/cost.py | 24 ++++++++++++++---------- tests/test_cost.py | 43 ++++++++++++++++++++++++++++++++++++++----- 2 files changed, 52 insertions(+), 15 deletions(-) diff --git a/src/iminuit/cost.py b/src/iminuit/cost.py index 95c75913..cb03fdaf 100644 --- a/src/iminuit/cost.py +++ b/src/iminuit/cost.py @@ -158,13 +158,13 @@ def __init__(self, val: ArrayLike, var: ArrayLike): @overload def __call__(self, val: ArrayLike) -> Tuple[NDArray, NDArray]: - ... + ... # pragma: no cover @overload def __call__( self, val: ArrayLike, var: ArrayLike ) -> Tuple[NDArray, NDArray, NDArray]: - ... + ... # pragma: no cover def __call__(self, val, var=None): """ @@ -181,11 +181,11 @@ def __call__(self, val, var=None): ------- (obs, pred) or (obs, pred, pred_var) """ + val = np.atleast_1d(val) s = self._scale - assert isinstance(val, np.ndarray) if var is None: return self._obs, val * s - assert isinstance(var, np.ndarray) + var = np.atleast_1d(var) return self._obs, val * s, var * s**2 @@ -590,7 +590,7 @@ def __call__(self, *args: float) -> float: @abc.abstractmethod def _call(self, args: Sequence[float]) -> float: - ... + ... # pragma: no cover class Constant(Cost): @@ -1272,12 +1272,16 @@ def __init__( t2 *= f**2 self._model_data.append((t1, t2)) args.append(f"x{i}") - else: - assert isinstance(t, Model) + elif isinstance(t, Model): par = describe(t)[1:] npar = len(par) self._model_data.append((t, npar)) args += [f"x{i}_{x}" for x in par] + else: + raise ValueError( + "model_or_template must be a collection of array-likes " + "and/or Model types" + ) if name is None: name = args @@ -1334,8 +1338,8 @@ def _pred(self, args: Sequence[float]) -> Tuple[NDArray, NDArray]: d = _normalize_model_output(d) if self._xe_shape is not None: d = d.reshape(self._xe_shape) - for i in range(self._ndim): - d = np.diff(d, axis=i) + for j in range(self._ndim): + d = np.diff(d, axis=j) # differences can come out negative due to round-off error in subtraction, # we set negative values to zero d[d < 0] = 0 @@ -1629,7 +1633,7 @@ def loss(self, loss: Union[str, LossFunction]): elif loss == "soft_l1": self._cost = _soft_l1_cost else: - raise ValueError("unknown loss type: " + loss) + raise ValueError(f"unknown loss {loss!r}") elif isinstance(loss, LossFunction): self._cost = lambda y, ye, ym: np.sum( loss(_z_squared(y, ye, ym)) # type:ignore diff --git a/tests/test_cost.py b/tests/test_cost.py index 2efb3845..95b1b2b9 100644 --- a/tests/test_cost.py +++ b/tests/test_cost.py @@ -733,15 +733,18 @@ def model(xy, a, b): def test_LeastSquares_bad_input(): - with pytest.raises(ValueError): + with pytest.raises(ValueError, match="shape mismatch"): LeastSquares([1, 2], [], [1], lambda x, a: 0) - with pytest.raises(ValueError): + with pytest.raises(ValueError, match="shape mismatch"): LeastSquares([1, 2], [3, 4, 5], [1], lambda x, a: 0) - with pytest.raises(ValueError): + with pytest.raises(ValueError, match="unknown loss"): LeastSquares([1], [1], [1], lambda x, a: 0, loss="foo") + with pytest.raises(ValueError, match="loss must be str or LossFunction"): + LeastSquares([1], [1], [1], lambda x, a: 0, loss=[1, 2, 3]) + def test_LeastSquares_mask(): c = LeastSquares([1, 2, 3], [3, np.nan, 4], [1, 1, 1], lambda x, a: x + a) @@ -1194,6 +1197,9 @@ def test_Template_bad_input(): with pytest.raises(ValueError, match="number of names"): Template([1], [1, 2], [[1]], name=("b", "s")) + with pytest.raises(ValueError, match="model_or_template"): + Template([1], [1, 2], [[1], None]) + @pytest.mark.skipif(not matplotlib_available, reason="matplotlib is needed") def test_Template_visualize(): @@ -1233,7 +1239,7 @@ def test_Template_pickle(): assert_equal(c.data, c2.data) -def test_Template_with_model_template_mix(): +def test_Template_with_model(): n = np.array([3, 2, 3]) xe = np.array([0, 1, 2, 3]) t = np.array([0.1, 0, 1]) @@ -1249,7 +1255,34 @@ def test_Template_with_model_template_mix(): assert m.valid -def test_Template_with_models(): +@pytest.mark.skipif(not scipy_available, reason="scipy.stats is needed") +def test_Template_with_model_2D(): + truth1 = (1.0, 0.1, 0.2, 0.3, 0.4, 0.5) + x1, y1 = mvnorm(*truth1[1:]).rvs(size=int(truth1[0] * 1000), random_state=1).T + truth2 = (1.0, 0.2, 0.1, 0.4, 0.3, 0.0) + x2, y2 = mvnorm(*truth2[1:]).rvs(size=int(truth2[0] * 1000), random_state=1).T + + x = np.append(x1, x2) + y = np.append(y1, y2) + w, xe, ye = np.histogram2d(x, y, bins=(3, 5)) + + def model(xy, n, mux, muy, sx, sy, rho): + return n * 1000 * mvnorm(mux, muy, sx, sy, rho).cdf(np.transpose(xy)) + + x3, y3 = mvnorm(*truth2[1:]).rvs(size=int(truth2[0] * 10000), random_state=2).T + template = np.histogram2d(x3, y3, bins=(xe, ye))[0] + + cost = Template(w, (xe, ye), (model, template)) + assert cost.ndata == np.prod(w.shape) + m = Minuit(cost, *truth1, 1) + m.limits["x0_n", "x0_sx", "x0_sy"] = (0, None) + m.limits["x0_rho"] = (-1, 1) + m.migrad() + assert m.valid + assert_allclose(m.values, truth1 + (1e3,), rtol=0.1) + + +def test_Template_with_only_models(): n = np.array([3, 2, 3]) xe = np.array([0, 1, 2, 3])