Nimare fix

numba/cpython/randomimpl.py

@@ -1782,6 +1782,13 @@ def multinomial_inner(n, pvals, out):
         fl = out.flat
         sz = out.size
         plen = len(pvals)
+
+        psum = 0.0
+        for i in range(plen):
+            psum += pvals[i]
+        if psum > 1 + 1e-10:
+            raise AssertionError("np.random.multinomial(): Sum of probabilities "
+                                 "cannot exceed 1.")
 
         for i in range(0, sz, plen):
             # Loop body: take a set of n experiments and fill up
@@ -1797,7 +1804,7 @@ def multinomial_inner(n, pvals, out):
             # distribution over the remaining number of experiments.
             for j in range(0, plen - 1):
                 p_j = pvals[j]
-                n_j = fl[i + j] = np.random.binomial(n_experiments, p_j / p_sum)
+                n_j = fl[i + j] = np.random.binomial(n_experiments, min(p_j / p_sum, 1))
                 n_experiments -= n_j
                 if n_experiments <= 0:
                     # Note the output was initialized to zero

numba/tests/test_random.py

@@ -1407,6 +1407,19 @@ def test_multinomial_3_tuple(self):
         for sample in res.reshape((-1, res.shape[-1])):
             self._check_sample(n, pvals, sample)
 
+    def test_multinomial_4_floating_point(self):
+        """
+        Test multinomial floating point
+        """
+        cfunc = jit(nopython=True)(numpy_multinomial2)
+        # The following values _just_ exceed one in sum
+        n, pvals = 2, np.array([0.4166666666666667, 0.5833333333333334, 0.0])
+        res = cfunc(n, pvals)
+        self.assertEqual(len(res), 3)
+        self.assertEqual(res[2], 0)
+        self.assertLessEqual(res[0], 2)
+        self.assertLessEqual(res[1], 2)
+
 
 class TestRandomDirichlet(BaseTest):
     alpha = np.array([1, 1, 1, 2], dtype=np.float64)