Skip to content

Commit dc441ee

Browse files
WarrenWeckessercharris
WarrenWeckesser
authored and
charris
committed
BUG: random: Create a legacy implementation of random.binomial.
Create a legacy implementation of the random_binomial method of RandomState that does not include the "short-circuit" check for n == 0 or p == 0. This ensures that the stream of variates is consistent with the behavior in 1.16. Closes numpygh-14522.
1 parent df93b2b commit dc441ee

7 files changed

Lines changed: 75 additions & 11 deletions

File tree

‎numpy/random/legacy_distributions.pxd‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -34,6 +34,8 @@ cdef extern from "legacy-distributions.h":
3434
double nonc) nogil
3535
double legacy_wald(aug_bitgen_t *aug_state, double mean, double scale) nogil
3636
double legacy_lognormal(aug_bitgen_t *aug_state, double mean, double sigma) nogil
37+
int64_t legacy_random_binomial(bitgen_t *bitgen_state, double p,
38+
int64_t n, binomial_t *binomial) nogil
3739
int64_t legacy_negative_binomial(aug_bitgen_t *aug_state, double n, double p) nogil
3840
int64_t legacy_random_hypergeometric(bitgen_t *bitgen_state, int64_t good, int64_t bad, int64_t sample) nogil
3941
int64_t legacy_random_logseries(bitgen_t *bitgen_state, double p) nogil

‎numpy/random/mtrand.pyx‎

Lines changed: 7 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -3086,7 +3086,9 @@ cdef class RandomState:
30863086
for i in range(cnt):
30873087
_dp = (<double*>np.PyArray_MultiIter_DATA(it, 1))[0]
30883088
_in = (<long*>np.PyArray_MultiIter_DATA(it, 2))[0]
3089-
(<long*>np.PyArray_MultiIter_DATA(it, 0))[0] = random_binomial(&self._bitgen, _dp, _in, &self._binomial)
3089+
(<long*>np.PyArray_MultiIter_DATA(it, 0))[0] = \
3090+
legacy_random_binomial(&self._bitgen, _dp, _in,
3091+
&self._binomial)
30903092

30913093
np.PyArray_MultiIter_NEXT(it)
30923094

@@ -3099,16 +3101,17 @@ cdef class RandomState:
30993101

31003102
if size is None:
31013103
with self.lock:
3102-
return random_binomial(&self._bitgen, _dp, _in, &self._binomial)
3104+
return <long>legacy_random_binomial(&self._bitgen, _dp, _in,
3105+
&self._binomial)
31033106

31043107
randoms = <np.ndarray>np.empty(size, int)
31053108
cnt = np.PyArray_SIZE(randoms)
31063109
randoms_data = <long *>np.PyArray_DATA(randoms)
31073110

31083111
with self.lock, nogil:
31093112
for i in range(cnt):
3110-
randoms_data[i] = random_binomial(&self._bitgen, _dp, _in,
3111-
&self._binomial)
3113+
randoms_data[i] = legacy_random_binomial(&self._bitgen, _dp, _in,
3114+
&self._binomial)
31123115

31133116
return randoms
31143117

‎numpy/random/src/distributions/distributions.c‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -901,8 +901,8 @@ RAND_INT_TYPE random_binomial_inversion(bitgen_t *bitgen_state, RAND_INT_TYPE n,
901901
return X;
902902
}
903903

904-
RAND_INT_TYPE random_binomial(bitgen_t *bitgen_state, double p, RAND_INT_TYPE n,
905-
binomial_t *binomial) {
904+
int64_t random_binomial(bitgen_t *bitgen_state, double p, int64_t n,
905+
binomial_t *binomial) {
906906
double q;
907907

908908
if ((n == 0LL) || (p == 0.0f))

‎numpy/random/src/distributions/distributions.h‎

Lines changed: 14 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -43,11 +43,11 @@
4343
typedef struct s_binomial_t {
4444
int has_binomial; /* !=0: following parameters initialized for binomial */
4545
double psave;
46-
int64_t nsave;
46+
RAND_INT_TYPE nsave;
4747
double r;
4848
double q;
4949
double fm;
50-
int64_t m;
50+
RAND_INT_TYPE m;
5151
double p1;
5252
double xm;
5353
double xl;
@@ -148,8 +148,18 @@ DECLDIR double random_triangular(bitgen_t *bitgen_state, double left, double mod
148148
DECLDIR RAND_INT_TYPE random_poisson(bitgen_t *bitgen_state, double lam);
149149
DECLDIR RAND_INT_TYPE random_negative_binomial(bitgen_t *bitgen_state, double n,
150150
double p);
151-
DECLDIR RAND_INT_TYPE random_binomial(bitgen_t *bitgen_state, double p, RAND_INT_TYPE n,
152-
binomial_t *binomial);
151+
152+
DECLDIR RAND_INT_TYPE random_binomial_btpe(bitgen_t *bitgen_state,
153+
RAND_INT_TYPE n,
154+
double p,
155+
binomial_t *binomial);
156+
DECLDIR RAND_INT_TYPE random_binomial_inversion(bitgen_t *bitgen_state,
157+
RAND_INT_TYPE n,
158+
double p,
159+
binomial_t *binomial);
160+
DECLDIR int64_t random_binomial(bitgen_t *bitgen_state, double p,
161+
int64_t n, binomial_t *binomial);
162+
153163
DECLDIR RAND_INT_TYPE random_logseries(bitgen_t *bitgen_state, double p);
154164
DECLDIR RAND_INT_TYPE random_geometric_search(bitgen_t *bitgen_state, double p);
155165
DECLDIR RAND_INT_TYPE random_geometric_inversion(bitgen_t *bitgen_state, double p);

‎numpy/random/src/legacy/legacy-distributions.c‎

Lines changed: 31 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -215,6 +215,37 @@ double legacy_exponential(aug_bitgen_t *aug_state, double scale) {
215215
}
216216

217217

218+
static RAND_INT_TYPE legacy_random_binomial_original(bitgen_t *bitgen_state,
219+
double p,
220+
RAND_INT_TYPE n,
221+
binomial_t *binomial) {
222+
double q;
223+
224+
if (p <= 0.5) {
225+
if (p * n <= 30.0) {
226+
return random_binomial_inversion(bitgen_state, n, p, binomial);
227+
} else {
228+
return random_binomial_btpe(bitgen_state, n, p, binomial);
229+
}
230+
} else {
231+
q = 1.0 - p;
232+
if (q * n <= 30.0) {
233+
return n - random_binomial_inversion(bitgen_state, n, q, binomial);
234+
} else {
235+
return n - random_binomial_btpe(bitgen_state, n, q, binomial);
236+
}
237+
}
238+
}
239+
240+
241+
int64_t legacy_random_binomial(bitgen_t *bitgen_state, double p,
242+
int64_t n, binomial_t *binomial) {
243+
return (int64_t) legacy_random_binomial_original(bitgen_state, p,
244+
(RAND_INT_TYPE) n,
245+
binomial);
246+
}
247+
248+
218249
static RAND_INT_TYPE random_hypergeometric_hyp(bitgen_t *bitgen_state,
219250
RAND_INT_TYPE good,
220251
RAND_INT_TYPE bad,

‎numpy/random/src/legacy/legacy-distributions.h‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -36,6 +36,8 @@ extern double legacy_f(aug_bitgen_t *aug_state, double dfnum, double dfden);
3636
extern double legacy_normal(aug_bitgen_t *aug_state, double loc, double scale);
3737
extern double legacy_standard_gamma(aug_bitgen_t *aug_state, double shape);
3838
extern double legacy_exponential(aug_bitgen_t *aug_state, double scale);
39+
extern int64_t legacy_random_binomial(bitgen_t *bitgen_state, double p,
40+
int64_t n, binomial_t *binomial);
3941
extern int64_t legacy_negative_binomial(aug_bitgen_t *aug_state, double n,
4042
double p);
4143
extern int64_t legacy_random_hypergeometric(bitgen_t *bitgen_state,

‎numpy/random/tests/test_randomstate_regression.py‎

Lines changed: 17 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -191,4 +191,20 @@ def test_randint_117(self):
191191
2588848963, 3684848379, 2340255427, 3638918503,
192192
1819583497, 2678185683], dtype='int64')
193193
actual = random.randint(2**32, size=10)
194-
assert_array_equal(actual, expected)
194+
assert_array_equal(actual, expected)
195+
196+
def test_p_zero_stream(self):
197+
# Regression test for gh-14522. Ensure that future versions
198+
# generate the same variates as version 1.16.
199+
np.random.seed(12345)
200+
assert_array_equal(random.binomial(1, [0, 0.25, 0.5, 0.75, 1]),
201+
[0, 0, 0, 1, 1])
202+
203+
def test_n_zero_stream(self):
204+
# Regression test for gh-14522. Ensure that future versions
205+
# generate the same variates as version 1.16.
206+
np.random.seed(8675309)
207+
expected = np.array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
208+
[3, 4, 2, 3, 3, 1, 5, 3, 1, 3]])
209+
assert_array_equal(random.binomial([[0], [10]], 0.25, size=(2, 10)),
210+
expected)

0 commit comments

Comments
 (0)