diff --git a/numpy/random/randomgen/common.pyx b/numpy/random/randomgen/common.pyx
index 1f7cd40cac76..ebd2e60d107d 100644
--- a/numpy/random/randomgen/common.pyx
+++ b/numpy/random/randomgen/common.pyx
@@ -337,6 +337,9 @@ cdef int check_constraint(double val, object name, constraint_type cons) except
     elif cons == CONS_BOUNDED_0_1:
         if not (val >= 0) or not (val <= 1):
             raise ValueError("{0} < 0 , {0} > 1 or {0} is NaN".format(name))
+    elif cons == CONS_BOUNDED_GT_0_1:
+        if not val >0 or not val <= 1:
+            raise ValueError("{0} <= 0 , {0} > 1 or {0} contains NaNs".format(name))
     elif cons == CONS_GT_1:
         if not (val > 1):
             raise ValueError("{0} <= 1 or {0} is NaN".format(name))
diff --git a/numpy/random/randomgen/src/legacy/distributions-boxmuller.c b/numpy/random/randomgen/src/legacy/distributions-boxmuller.c
index e96d3f19a0cb..5d3ba27f8218 100644
--- a/numpy/random/randomgen/src/legacy/distributions-boxmuller.c
+++ b/numpy/random/randomgen/src/legacy/distributions-boxmuller.c
@@ -12,7 +12,6 @@ double legacy_gauss(aug_brng_t *aug_state) {
     return temp;
   } else {
     double f, x1, x2, r2;
-    double f, x1, x2, r2;
 
     do {
       x1 = 2.0 * legacy_double(aug_state) - 1.0;
diff --git a/numpy/random/randomgen/tests/test_against_numpy.py b/numpy/random/randomgen/tests/test_against_numpy.py
index b0fcdfb710f1..431c7bd85d24 100644
--- a/numpy/random/randomgen/tests/test_against_numpy.py
+++ b/numpy/random/randomgen/tests/test_against_numpy.py
@@ -570,9 +570,9 @@ def test_multivariate_normal(self):
         self._is_state_common_legacy()
 
 
-funcs = [randomgen.generator.zipf,
-         randomgen.generator.logseries,
-         randomgen.generator.poisson]
+funcs = [generator.zipf,
+         generator.logseries,
+         generator.poisson]
 ids = [f.__name__ for f in funcs]
 
 
diff --git a/numpy/random/randomgen/tests/test_direct.py b/numpy/random/randomgen/tests/test_direct.py
index 6e856de414cb..22fdcd865944 100644
--- a/numpy/random/randomgen/tests/test_direct.py
+++ b/numpy/random/randomgen/tests/test_direct.py
@@ -233,14 +233,14 @@ def test_seed_out_of_range_array(self):
 
     def test_repr(self):
         rs = RandomGenerator(self.brng(*self.data1['seed']))
-        assert 'RandomGenerator' in rs.__repr__()
-        assert str(hex(id(rs)))[2:].upper() in rs.__repr__()
+        assert 'RandomGenerator' in repr(rs)
+        assert '{:#x}'.format(id(rs)).upper().replace('X', 'x') in repr(rs)
 
     def test_str(self):
         rs = RandomGenerator(self.brng(*self.data1['seed']))
         assert 'RandomGenerator' in str(rs)
         assert str(self.brng.__name__) in str(rs)
-        assert str(hex(id(rs)))[2:].upper() not in str(rs)
+        assert '{:#x}'.format(id(rs)).upper().replace('X', 'x') not in str(rs)
 
     def test_generator(self):
         brng = self.brng(*self.data1['seed'])
diff --git a/numpy/random/randomgen/tests/test_generator_mt19937.py b/numpy/random/randomgen/tests/test_generator_mt19937.py
index cad3ef4d66cd..10dd8733a8b9 100644
--- a/numpy/random/randomgen/tests/test_generator_mt19937.py
+++ b/numpy/random/randomgen/tests/test_generator_mt19937.py
@@ -53,8 +53,7 @@ def test_n_zero(self):
         # This test addresses issue #3480.
         zeros = np.zeros(2, dtype='int')
         for p in [0, .5, 1]:
-            val = random.binomial(0, p)
-            assert val == 0
+            assert_(random.binomial(0, p) == 0)
             assert_array_equal(random.binomial(zeros, p), zeros)
 
     def test_p_is_nan(self):
@@ -96,40 +95,40 @@ def test_invalid_prob(self):
 class TestSetState(object):
     def setup(self):
         self.seed = 1234567890
-        self.brng = RandomGenerator(MT19937(self.seed))
-        self.state = self.brng.state
+        self.rg = RandomGenerator(MT19937(self.seed))
+        self.state = self.rg.state
         self.legacy_state = (self.state['brng'],
                              self.state['state']['key'],
                              self.state['state']['pos'])
 
     def test_basic(self):
-        old = self.brng.tomaxint(16)
-        self.brng.state = self.state
-        new = self.brng.tomaxint(16)
+        old = self.rg.tomaxint(16)
+        self.rg.state = self.state
+        new = self.rg.tomaxint(16)
         assert_(np.all(old == new))
 
     def test_gaussian_reset(self):
         # Make sure the cached every-other-Gaussian is reset.
-        old = self.brng.standard_normal(size=3)
-        self.brng.state = self.state
-        new = self.brng.standard_normal(size=3)
+        old = self.rg.standard_normal(size=3)
+        self.rg.state = self.state
+        new = self.rg.standard_normal(size=3)
         assert_(np.all(old == new))
 
     def test_gaussian_reset_in_media_res(self):
         # When the state is saved with a cached Gaussian, make sure the
         # cached Gaussian is restored.
 
-        self.brng.standard_normal()
-        state = self.brng.state
-        old = self.brng.standard_normal(size=3)
-        self.brng.state = state
-        new = self.brng.standard_normal(size=3)
+        self.rg.standard_normal()
+        state = self.rg.state
+        old = self.rg.standard_normal(size=3)
+        self.rg.state = state
+        new = self.rg.standard_normal(size=3)
         assert_(np.all(old == new))
 
     def test_negative_binomial(self):
         # Ensure that the negative binomial results take floating point
         # arguments without truncation.
-        self.brng.negative_binomial(0.5, 0.5)
+        self.rg.negative_binomial(0.5, 0.5)
 
 
 class TestRandint(object):
@@ -554,8 +553,7 @@ def test_choice_uniform_noreplace(self):
 
     def test_choice_nonuniform_noreplace(self):
         random.seed(self.seed)
-        actual = random.choice(4, 3, replace=False,
-                               p=[0.1, 0.3, 0.5, 0.1])
+        actual = random.choice(4, 3, replace=False, p=[0.1, 0.3, 0.5, 0.1])
         desired = np.array([2, 3, 1])
         assert_array_equal(actual, desired)
 
@@ -614,11 +612,11 @@ def test_choice_return_shape(self):
         # Check multi dimensional array
         s = (2, 3)
         p = [0.1, 0.1, 0.1, 0.1, 0.4, 0.2]
-        assert_(random.choice(6, s, replace=True).shape, s)
-        assert_(random.choice(6, s, replace=False).shape, s)
-        assert_(random.choice(6, s, replace=True, p=p).shape, s)
-        assert_(random.choice(6, s, replace=False, p=p).shape, s)
-        assert_(random.choice(np.arange(6), s, replace=True).shape, s)
+        assert_equal(random.choice(6, s, replace=True).shape, s)
+        assert_equal(random.choice(6, s, replace=False).shape, s)
+        assert_equal(random.choice(6, s, replace=True, p=p).shape, s)
+        assert_equal(random.choice(6, s, replace=False, p=p).shape, s)
+        assert_equal(random.choice(np.arange(6), s, replace=True).shape, s)
 
         # Check zero-size
         assert_equal(random.randint(0, 0, size=(3, 0, 4)).shape, (3, 0, 4))
@@ -711,6 +709,11 @@ def test_binomial(self):
                             [46, 45]])
         assert_array_equal(actual, desired)
 
+        random.seed(self.seed)
+        actual = random.binomial(100.123, .456)
+        desired = 37
+        assert_array_equal(actual, desired)
+
     def test_chisquare(self):
         random.seed(self.seed)
         actual = random.chisquare(50, size=(3, 2))
@@ -795,6 +798,16 @@ def test_geometric(self):
                             [5, 12]])
         assert_array_equal(actual, desired)
 
+    def test_geometric_exceptions(self):
+        assert_raises(ValueError, random.geometric, 1.1)
+        assert_raises(ValueError, random.geometric, [1.1] * 10)
+        assert_raises(ValueError, random.geometric, -0.1)
+        assert_raises(ValueError, random.geometric, [-0.1] * 10)
+        with suppress_warnings() as sup:
+            sup.record(RuntimeWarning)
+            assert_raises(ValueError, random.geometric, np.nan)
+            assert_raises(ValueError, random.geometric, [np.nan] * 10)
+
     def test_gumbel(self):
         random.seed(self.seed)
         actual = random.gumbel(loc=.123456789, scale=2.0, size=(3, 2))
@@ -873,6 +886,12 @@ def test_logseries(self):
                             [3, 6]])
         assert_array_equal(actual, desired)
 
+    def test_logseries_exceptions(self):
+        with suppress_warnings() as sup:
+            sup.record(RuntimeWarning)
+            assert_raises(ValueError, random.logseries, np.nan)
+            assert_raises(ValueError, random.logseries, [np.nan] * 10)
+
     def test_multinomial(self):
         random.seed(self.seed)
         actual = random.multinomial(20, [1 / 6.] * 6, size=(3, 2))
@@ -943,6 +962,13 @@ def test_negative_binomial(self):
                             [723, 751]])
         assert_array_equal(actual, desired)
 
+    def test_negative_binomial_exceptions(self):
+        with suppress_warnings() as sup:
+            sup.record(RuntimeWarning)
+            assert_raises(ValueError, random.negative_binomial, 100, np.nan)
+            assert_raises(ValueError, random.negative_binomial, 100,
+                          [np.nan] * 10)
+
     def test_noncentral_chisquare(self):
         random.seed(self.seed)
         actual = random.noncentral_chisquare(df=5, nonc=5, size=(3, 2))
@@ -973,6 +999,11 @@ def test_noncentral_f(self):
                             [1.16362730891403, 2.54104276581491]])
         assert_array_almost_equal(actual, desired, decimal=14)
 
+    def test_noncentral_f_nan(self):
+        random.seed(self.seed)
+        actual = random.noncentral_f(dfnum=5, dfden=2, nonc=np.nan)
+        assert np.isnan(actual)
+
     def test_normal(self):
         random.seed(self.seed)
         actual = random.normal(loc=.123456789, scale=2.0, size=(3, 2))
@@ -993,7 +1024,7 @@ def test_pareto(self):
                             [1.1281132447159091e+01, 3.1895968171107006e+08]])
         # For some reason on 32-bit x86 Ubuntu 12.10 the [1, 0] entry in this
         # matrix differs by 24 nulps. Discussion:
-        #   http://mail.scipy.org/pipermail/numpy-discussion/2012-September/063801.html
+        #   https://mail.python.org/pipermail/numpy-discussion/2012-September/063801.html
         # Consensus is that this is probably some gcc quirk that affects
         # rounding but not in any important way, so we just use a looser
         # tolerance on this test:
@@ -1014,6 +1045,10 @@ def test_poisson_exceptions(self):
         assert_raises(ValueError, random.poisson, [lamneg] * 10)
         assert_raises(ValueError, random.poisson, lambig)
         assert_raises(ValueError, random.poisson, [lambig] * 10)
+        with suppress_warnings() as sup:
+            sup.record(RuntimeWarning)
+            assert_raises(ValueError, random.poisson, np.nan)
+            assert_raises(ValueError, random.poisson, [np.nan] * 10)
 
     def test_power(self):
         random.seed(self.seed)
@@ -1168,8 +1203,8 @@ def __float__(self):
                 raise TypeError
 
         throwing_float = np.array(1.0).view(ThrowingFloat)
-        assert_raises(TypeError, random.uniform,
-                      throwing_float, throwing_float)
+        assert_raises(TypeError, random.uniform, throwing_float,
+                      throwing_float)
 
         class ThrowingInteger(np.ndarray):
             def __int__(self):
@@ -1189,9 +1224,14 @@ def test_vonmises(self):
     def test_vonmises_small(self):
         # check infinite loop, gh-4720
         random.seed(self.seed)
-        r = random.vonmises(mu=0., kappa=1.1e-8, size=10 ** 6)
+        r = random.vonmises(mu=0., kappa=1.1e-8, size=10**6)
         assert_(np.isfinite(r).all())
 
+    def test_vonmises_nan(self):
+        random.seed(self.seed)
+        r = random.vonmises(mu=0., kappa=np.nan)
+        assert_(np.isnan(r))
+
     def test_wald(self):
         random.seed(self.seed)
         actual = random.wald(mean=1.23, scale=1.54, size=(3, 2))
@@ -1372,8 +1412,9 @@ def test_noncentral_f(self):
 
         self.set_seed()
         actual = nonc_f(dfnum * 3, dfden, nonc)
-        mt_nonc_f = random.noncentral_f
         assert_array_almost_equal(actual, desired, decimal=14)
+        assert np.all(np.isnan(nonc_f(dfnum, dfden, [np.nan] * 3)))
+
         assert_raises(ValueError, nonc_f, bad_dfnum * 3, dfden, nonc)
         assert_raises(ValueError, nonc_f, dfnum * 3, bad_dfden, nonc)
         assert_raises(ValueError, nonc_f, dfnum * 3, dfden, bad_nonc)
@@ -1391,9 +1432,12 @@ def test_noncentral_f(self):
         assert_raises(ValueError, nonc_f, bad_dfnum, dfden, nonc * 3)
         assert_raises(ValueError, nonc_f, dfnum, bad_dfden, nonc * 3)
         assert_raises(ValueError, nonc_f, dfnum, dfden, bad_nonc * 3)
-        assert_raises(ValueError, mt_nonc_f, bad_dfnum, dfden, nonc * 3)
-        assert_raises(ValueError, mt_nonc_f, dfnum, bad_dfden, nonc * 3)
-        assert_raises(ValueError, mt_nonc_f, dfnum, dfden, bad_nonc * 3)
+
+    def test_noncentral_f_small_df(self):
+        self.set_seed()
+        desired = np.array([21.57878070681719,  1.17110217503908])
+        actual = random.noncentral_f(0.9, 0.9, 2, size=2)
+        assert_array_almost_equal(actual, desired, decimal=14)
 
     def test_chisquare(self):
         df = [1]
@@ -1420,20 +1464,15 @@ def test_noncentral_chisquare(self):
 
         self.set_seed()
         actual = nonc_chi(df * 3, nonc)
-        mt_nonc_chi2 = random.noncentral_chisquare
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, nonc_chi, bad_df * 3, nonc)
         assert_raises(ValueError, nonc_chi, df * 3, bad_nonc)
-        assert_raises(ValueError, mt_nonc_chi2, bad_df * 3, nonc)
-        assert_raises(ValueError, mt_nonc_chi2, df * 3, bad_nonc)
 
         self.set_seed()
         actual = nonc_chi(df, nonc * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, nonc_chi, bad_df, nonc * 3)
         assert_raises(ValueError, nonc_chi, df, bad_nonc * 3)
-        assert_raises(ValueError, mt_nonc_chi2, bad_df, nonc * 3)
-        assert_raises(ValueError, mt_nonc_chi2, df, bad_nonc * 3)
 
     def test_standard_t(self):
         df = [1]
@@ -1645,24 +1684,24 @@ def test_triangular(self):
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, triangular, bad_left_one * 3, mode, right)
         assert_raises(ValueError, triangular, left * 3, bad_mode_one, right)
-        assert_raises(ValueError, triangular,
-                      bad_left_two * 3, bad_mode_two, right)
+        assert_raises(ValueError, triangular, bad_left_two * 3, bad_mode_two,
+                      right)
 
         self.set_seed()
         actual = triangular(left, mode * 3, right)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, triangular, bad_left_one, mode * 3, right)
         assert_raises(ValueError, triangular, left, bad_mode_one * 3, right)
-        assert_raises(ValueError, triangular, bad_left_two,
-                      bad_mode_two * 3, right)
+        assert_raises(ValueError, triangular, bad_left_two, bad_mode_two * 3,
+                      right)
 
         self.set_seed()
         actual = triangular(left, mode, right * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, triangular, bad_left_one, mode, right * 3)
         assert_raises(ValueError, triangular, left, bad_mode_one, right * 3)
-        assert_raises(ValueError, triangular, bad_left_two,
-                      bad_mode_two, right * 3)
+        assert_raises(ValueError, triangular, bad_left_two, bad_mode_two,
+                      right * 3)
 
         assert_raises(ValueError, triangular, 10., 0., 20.)
         assert_raises(ValueError, triangular, 10., 25., 20.)
@@ -1739,6 +1778,9 @@ def test_zipf(self):
         actual = zipf(a * 3)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, zipf, bad_a * 3)
+        with np.errstate(invalid='ignore'):
+            assert_raises(ValueError, zipf, np.nan)
+            assert_raises(ValueError, zipf, [0, 0, np.nan])
 
     def test_geometric(self):
         p = [0.5]
@@ -1809,7 +1851,6 @@ def test_logseries(self):
 
 class TestThread(object):
     # make sure each state produces the same sequence even in threads
-
     def setup(self):
         self.seeds = range(4)
 
diff --git a/numpy/random/randomgen/tests/test_randomstate.py b/numpy/random/randomgen/tests/test_randomstate.py
index 167d1b0aa545..6adde6c2ab33 100644
--- a/numpy/random/randomgen/tests/test_randomstate.py
+++ b/numpy/random/randomgen/tests/test_randomstate.py
@@ -109,7 +109,9 @@ def test_size(self):
 
         assert_raises(TypeError, random.multinomial, 1, p,
                       float(1))
-        assert_raises(ValueError, random.multinomial, 1, [1.1, .1])
+
+    def test_invalid_prob(self):
+        assert_raises(ValueError, random.multinomial, 100, [1.1, 0.2])
 
 
 class TestSetState(object):
@@ -352,6 +354,12 @@ def test_rand(self):
                             [0.4575674820298663, 0.7781880808593471]])
         assert_array_almost_equal(actual, desired, decimal=15)
 
+    def test_rand_singleton(self):
+        random.seed(self.seed)
+        actual = random.rand()
+        desired = 0.61879477158567997
+        assert_array_almost_equal(actual, desired, decimal=15)
+
     def test_randn(self):
         random.seed(self.seed)
         actual = random.randn(3, 2)
@@ -360,6 +368,10 @@ def test_randn(self):
                            [2.031033998682787, 2.17032494605655257]])
         assert_array_almost_equal(actual, desired, decimal=15)
 
+        random.seed(self.seed)
+        actual = random.randn()
+        assert_array_almost_equal(actual, desired[0, 0], decimal=15)
+
     def test_randint(self):
         random.seed(self.seed)
         actual = random.randint(-99, 99, size=(3, 2))
@@ -443,11 +455,9 @@ def test_random_sample(self):
                             [0.4575674820298663, 0.7781880808593471]])
         assert_array_almost_equal(actual, desired, decimal=15)
 
-    def test_rand_singleton(self):
         random.seed(self.seed)
-        actual = random.rand()
-        desired = np.array(0.61879477158567997)
-        assert_array_almost_equal(actual, desired, decimal=15)
+        actual = random.random_sample()
+        assert_array_almost_equal(actual, desired[0, 0], decimal=15)
 
     def test_choice_uniform_replace(self):
         random.seed(self.seed)
@@ -582,18 +592,6 @@ def test_shuffle(self):
             desired = conv([0, 1, 9, 6, 2, 4, 5, 8, 7, 3])
             assert_array_equal(actual, desired)
 
-    def test_permutation(self):
-        random.seed(self.seed)
-        alist = [1, 2, 3, 4, 5, 6, 7, 8, 9, 0]
-        actual = random.permutation(alist)
-        desired = [0, 1, 9, 6, 2, 4, 5, 8, 7, 3]
-        assert_array_equal(actual, desired)
-
-        random.seed(self.seed)
-        arr_2d = np.atleast_2d([1, 2, 3, 4, 5, 6, 7, 8, 9, 0]).T
-        actual = random.permutation(arr_2d)
-        assert_array_equal(actual, np.atleast_2d(desired).T)
-
     def test_shuffle_masked(self):
         # gh-3263
         a = np.ma.masked_values(np.reshape(range(20), (5, 4)) % 3 - 1, -1)
@@ -608,6 +606,18 @@ def test_shuffle_masked(self):
             assert_equal(
                 sorted(b.data[~b.mask]), sorted(b_orig.data[~b_orig.mask]))
 
+    def test_permutation(self):
+        random.seed(self.seed)
+        alist = [1, 2, 3, 4, 5, 6, 7, 8, 9, 0]
+        actual = random.permutation(alist)
+        desired = [0, 1, 9, 6, 2, 4, 5, 8, 7, 3]
+        assert_array_equal(actual, desired)
+
+        random.seed(self.seed)
+        arr_2d = np.atleast_2d([1, 2, 3, 4, 5, 6, 7, 8, 9, 0]).T
+        actual = random.permutation(arr_2d)
+        assert_array_equal(actual, np.atleast_2d(desired).T)
+
     def test_beta(self):
         random.seed(self.seed)
         actual = random.beta(.1, .9, size=(3, 2))
@@ -649,6 +659,8 @@ def test_dirichlet(self):
                             [[0.59266909280647828, 0.40733090719352177],
                              [0.56974431743975207, 0.43025568256024799]]])
         assert_array_almost_equal(actual, desired, decimal=15)
+        bad_alpha = np.array([5.4e-01, -1.0e-16])
+        assert_raises(ValueError, random.dirichlet, bad_alpha)
 
         random.seed(self.seed)
         alpha = np.array([51.72840233779265162, 39.74494232180943953])
@@ -712,6 +724,16 @@ def test_geometric(self):
                             [5, 12]])
         assert_array_equal(actual, desired)
 
+    def test_geometric_exceptions(self):
+        assert_raises(ValueError, random.geometric, 1.1)
+        assert_raises(ValueError, random.geometric, [1.1] * 10)
+        assert_raises(ValueError, random.geometric, -0.1)
+        assert_raises(ValueError, random.geometric, [-0.1] * 10)
+        with suppress_warnings() as sup:
+            sup.record(RuntimeWarning)
+            assert_raises(ValueError, random.geometric, np.nan)
+            assert_raises(ValueError, random.geometric, [np.nan] * 10)
+
     def test_gumbel(self):
         random.seed(self.seed)
         actual = random.gumbel(loc=.123456789, scale=2.0, size=(3, 2))
@@ -790,9 +812,15 @@ def test_logseries(self):
                             [3, 6]])
         assert_array_equal(actual, desired)
 
+    def test_logseries_exceptions(self):
+        with suppress_warnings() as sup:
+            sup.record(RuntimeWarning)
+            assert_raises(ValueError, random.logseries, np.nan)
+            assert_raises(ValueError, random.logseries, [np.nan] * 10)
+
     def test_multinomial(self):
         random.seed(self.seed)
-        actual = random.multinomial(20, [1/6.]*6, size=(3, 2))
+        actual = random.multinomial(20, [1 / 6.] * 6, size=(3, 2))
         desired = np.array([[[4, 3, 5, 4, 2, 2],
                              [5, 2, 8, 2, 2, 1]],
                             [[3, 4, 3, 6, 0, 4],
@@ -843,8 +871,8 @@ def test_multivariate_normal(self):
 
         mu = np.zeros(2)
         cov = np.eye(2)
-        assert_raises(ValueError, random.multivariate_normal, mean,
-                      cov, check_valid='other')
+        assert_raises(ValueError, random.multivariate_normal, mean, cov,
+                      check_valid='other')
         assert_raises(ValueError, random.multivariate_normal,
                       np.zeros((2, 1, 1)), cov)
         assert_raises(ValueError, random.multivariate_normal,
@@ -860,6 +888,13 @@ def test_negative_binomial(self):
                             [779, 647]])
         assert_array_equal(actual, desired)
 
+    def test_negative_binomial_exceptions(self):
+        with suppress_warnings() as sup:
+            sup.record(RuntimeWarning)
+            assert_raises(ValueError, random.negative_binomial, 100, np.nan)
+            assert_raises(ValueError, random.negative_binomial, 100,
+                          [np.nan] * 10)
+
     def test_noncentral_chisquare(self):
         random.seed(self.seed)
         actual = random.noncentral_chisquare(df=5, nonc=5, size=(3, 2))
@@ -890,6 +925,11 @@ def test_noncentral_f(self):
                             [0.43741599463544162, 1.1774208752428319]])
         assert_array_almost_equal(actual, desired, decimal=14)
 
+    def test_noncentral_f_nan(self):
+        random.seed(self.seed)
+        actual = random.noncentral_f(dfnum=5, dfden=2, nonc=np.nan)
+        assert np.isnan(actual)
+
     def test_normal(self):
         random.seed(self.seed)
         actual = random.normal(loc=.123456789, scale=2.0, size=(3, 2))
@@ -929,9 +969,13 @@ def test_poisson_exceptions(self):
         lambig = np.iinfo('l').max
         lamneg = -1
         assert_raises(ValueError, random.poisson, lamneg)
-        assert_raises(ValueError, random.poisson, [lamneg]*10)
+        assert_raises(ValueError, random.poisson, [lamneg] * 10)
         assert_raises(ValueError, random.poisson, lambig)
-        assert_raises(ValueError, random.poisson, [lambig]*10)
+        assert_raises(ValueError, random.poisson, [lambig] * 10)
+        with suppress_warnings() as sup:
+            sup.record(RuntimeWarning)
+            assert_raises(ValueError, random.poisson, np.nan)
+            assert_raises(ValueError, random.poisson, [np.nan] * 10)
 
     def test_power(self):
         random.seed(self.seed)
@@ -1026,8 +1070,8 @@ def test_uniform_range_bounds(self):
 
         func = random.uniform
         assert_raises(OverflowError, func, -np.inf, 0)
-        assert_raises(OverflowError, func,  0,      np.inf)
-        assert_raises(OverflowError, func,  fmin,   fmax)
+        assert_raises(OverflowError, func, 0, np.inf)
+        assert_raises(OverflowError, func, fmin, fmax)
         assert_raises(OverflowError, func, [-np.inf], [0])
         assert_raises(OverflowError, func, [0], [np.inf])
 
@@ -1070,7 +1114,12 @@ def test_vonmises_small(self):
         # check infinite loop, gh-4720
         random.seed(self.seed)
         r = random.vonmises(mu=0., kappa=1.1e-8, size=10**6)
-        np.testing.assert_(np.isfinite(r).all())
+        assert_(np.isfinite(r).all())
+
+    def test_vonmises_nan(self):
+        random.seed(self.seed)
+        r = random.vonmises(mu=0., kappa=np.nan)
+        assert_(np.isnan(r))
 
     def test_wald(self):
         random.seed(self.seed)
@@ -1108,12 +1157,9 @@ class TestBroadcast(object):
     def setup(self):
         self.seed = 123456789
 
-    def setSeed(self):
+    def set_seed(self):
         random.seed(self.seed)
 
-    # TODO: Include test for randint once it can broadcast
-    # Can steal the test written in PR #6938
-
     def test_uniform(self):
         low = [0]
         high = [1]
@@ -1122,11 +1168,11 @@ def test_uniform(self):
                             0.53413660089041659,
                             0.50955303552646702])
 
-        self.setSeed()
+        self.set_seed()
         actual = uniform(low * 3, high)
         assert_array_almost_equal(actual, desired, decimal=14)
 
-        self.setSeed()
+        self.set_seed()
         actual = uniform(low, high * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
 
@@ -1139,12 +1185,12 @@ def test_normal(self):
                             2.1283977976520019,
                             1.8417114045748335])
 
-        self.setSeed()
+        self.set_seed()
         actual = normal(loc * 3, scale)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, normal, loc * 3, bad_scale)
 
-        self.setSeed()
+        self.set_seed()
         actual = normal(loc, scale * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, normal, loc, bad_scale * 3)
@@ -1159,13 +1205,13 @@ def test_beta(self):
                             0.075230336409423643,
                             0.24976865978980844])
 
-        self.setSeed()
+        self.set_seed()
         actual = beta(a * 3, b)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, beta, bad_a * 3, b)
         assert_raises(ValueError, beta, a * 3, bad_b)
 
-        self.setSeed()
+        self.set_seed()
         actual = beta(a, b * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, beta, bad_a, b * 3)
@@ -1179,7 +1225,7 @@ def test_exponential(self):
                             0.76386282278691653,
                             0.71243813125891797])
 
-        self.setSeed()
+        self.set_seed()
         actual = exponential(scale * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, exponential, bad_scale * 3)
@@ -1192,7 +1238,7 @@ def test_standard_gamma(self):
                             0.76386282278691653,
                             0.71243813125891797])
 
-        self.setSeed()
+        self.set_seed()
         actual = std_gamma(shape * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, std_gamma, bad_shape * 3)
@@ -1207,13 +1253,13 @@ def test_gamma(self):
                             1.5277256455738331,
                             1.4248762625178359])
 
-        self.setSeed()
+        self.set_seed()
         actual = gamma(shape * 3, scale)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, gamma, bad_shape * 3, scale)
         assert_raises(ValueError, gamma, shape * 3, bad_scale)
 
-        self.setSeed()
+        self.set_seed()
         actual = gamma(shape, scale * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, gamma, bad_shape, scale * 3)
@@ -1229,13 +1275,13 @@ def test_f(self):
                             0.86768719635363512,
                             2.7251095168386801])
 
-        self.setSeed()
+        self.set_seed()
         actual = f(dfnum * 3, dfden)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, f, bad_dfnum * 3, dfden)
         assert_raises(ValueError, f, dfnum * 3, bad_dfden)
 
-        self.setSeed()
+        self.set_seed()
         actual = f(dfnum, dfden * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, f, bad_dfnum, dfden * 3)
@@ -1253,21 +1299,23 @@ def test_noncentral_f(self):
                             13.025456344595602,
                             8.8018098359100545])
 
-        self.setSeed()
+        self.set_seed()
         actual = nonc_f(dfnum * 3, dfden, nonc)
         assert_array_almost_equal(actual, desired, decimal=14)
+        assert np.all(np.isnan(nonc_f(dfnum, dfden, [np.nan] * 3)))
+
         assert_raises(ValueError, nonc_f, bad_dfnum * 3, dfden, nonc)
         assert_raises(ValueError, nonc_f, dfnum * 3, bad_dfden, nonc)
         assert_raises(ValueError, nonc_f, dfnum * 3, dfden, bad_nonc)
 
-        self.setSeed()
+        self.set_seed()
         actual = nonc_f(dfnum, dfden * 3, nonc)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, nonc_f, bad_dfnum, dfden * 3, nonc)
         assert_raises(ValueError, nonc_f, dfnum, bad_dfden * 3, nonc)
         assert_raises(ValueError, nonc_f, dfnum, dfden * 3, bad_nonc)
 
-        self.setSeed()
+        self.set_seed()
         actual = nonc_f(dfnum, dfden, nonc * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, nonc_f, bad_dfnum, dfden, nonc * 3)
@@ -1275,7 +1323,7 @@ def test_noncentral_f(self):
         assert_raises(ValueError, nonc_f, dfnum, dfden, bad_nonc * 3)
 
     def test_noncentral_f_small_df(self):
-        self.setSeed()
+        self.set_seed()
         desired = np.array([6.869638627492048, 0.785880199263955])
         actual = random.noncentral_f(0.9, 0.9, 2, size=2)
         assert_array_almost_equal(actual, desired, decimal=14)
@@ -1288,7 +1336,7 @@ def test_chisquare(self):
                             0.51947702108840776,
                             0.1320969254923558])
 
-        self.setSeed()
+        self.set_seed()
         actual = chisquare(df * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, chisquare, bad_df * 3)
@@ -1303,13 +1351,13 @@ def test_noncentral_chisquare(self):
                             4.5804135049718742,
                             6.0872302432834564])
 
-        self.setSeed()
+        self.set_seed()
         actual = nonc_chi(df * 3, nonc)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, nonc_chi, bad_df * 3, nonc)
         assert_raises(ValueError, nonc_chi, df * 3, bad_nonc)
 
-        self.setSeed()
+        self.set_seed()
         actual = nonc_chi(df, nonc * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, nonc_chi, bad_df, nonc * 3)
@@ -1323,10 +1371,11 @@ def test_standard_t(self):
                             5.8560725167361607,
                             1.0274791436474273])
 
-        self.setSeed()
+        self.set_seed()
         actual = t(df * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, t, bad_df * 3)
+        assert_raises(ValueError, random.standard_t, bad_df * 3)
 
     def test_vonmises(self):
         mu = [2]
@@ -1337,12 +1386,12 @@ def test_vonmises(self):
                             -2.7064099483995943,
                             -1.8672476700665914])
 
-        self.setSeed()
+        self.set_seed()
         actual = vonmises(mu * 3, kappa)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, vonmises, mu * 3, bad_kappa)
 
-        self.setSeed()
+        self.set_seed()
         actual = vonmises(mu, kappa * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, vonmises, mu, bad_kappa * 3)
@@ -1355,10 +1404,11 @@ def test_pareto(self):
                             1.1465519762044529,
                             1.0389564467453547])
 
-        self.setSeed()
+        self.set_seed()
         actual = pareto(a * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, pareto, bad_a * 3)
+        assert_raises(ValueError, random.pareto, bad_a * 3)
 
     def test_weibull(self):
         a = [1]
@@ -1368,10 +1418,11 @@ def test_weibull(self):
                             0.76386282278691653,
                             0.71243813125891797])
 
-        self.setSeed()
+        self.set_seed()
         actual = weibull(a * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, weibull, bad_a * 3)
+        assert_raises(ValueError, random.weibull, bad_a * 3)
 
     def test_power(self):
         a = [1]
@@ -1381,10 +1432,11 @@ def test_power(self):
                             0.53413660089041659,
                             0.50955303552646702])
 
-        self.setSeed()
+        self.set_seed()
         actual = power(a * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, power, bad_a * 3)
+        assert_raises(ValueError, random.power, bad_a * 3)
 
     def test_laplace(self):
         loc = [0]
@@ -1395,12 +1447,12 @@ def test_laplace(self):
                             0.070715642226971326,
                             0.019290950698972624])
 
-        self.setSeed()
+        self.set_seed()
         actual = laplace(loc * 3, scale)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, laplace, loc * 3, bad_scale)
 
-        self.setSeed()
+        self.set_seed()
         actual = laplace(loc, scale * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, laplace, loc, bad_scale * 3)
@@ -1414,12 +1466,12 @@ def test_gumbel(self):
                             0.26936705726291116,
                             0.33906220393037939])
 
-        self.setSeed()
+        self.set_seed()
         actual = gumbel(loc * 3, scale)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, gumbel, loc * 3, bad_scale)
 
-        self.setSeed()
+        self.set_seed()
         actual = gumbel(loc, scale * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, gumbel, loc, bad_scale * 3)
@@ -1433,15 +1485,16 @@ def test_logistic(self):
                             0.13675915696285773,
                             0.038216792802833396])
 
-        self.setSeed()
+        self.set_seed()
         actual = logistic(loc * 3, scale)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, logistic, loc * 3, bad_scale)
 
-        self.setSeed()
+        self.set_seed()
         actual = logistic(loc, scale * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, logistic, loc, bad_scale * 3)
+        assert_equal(random.logistic(1.0, 0.0), 1.0)
 
     def test_lognormal(self):
         mean = [0]
@@ -1452,15 +1505,17 @@ def test_lognormal(self):
                             8.4013952870126261,
                             6.3073234116578671])
 
-        self.setSeed()
+        self.set_seed()
         actual = lognormal(mean * 3, sigma)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, lognormal, mean * 3, bad_sigma)
+        assert_raises(ValueError, random.lognormal, mean * 3, bad_sigma)
 
-        self.setSeed()
+        self.set_seed()
         actual = lognormal(mean, sigma * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, lognormal, mean, bad_sigma * 3)
+        assert_raises(ValueError, random.lognormal, mean, bad_sigma * 3)
 
     def test_rayleigh(self):
         scale = [1]
@@ -1470,7 +1525,7 @@ def test_rayleigh(self):
                             1.2360119924878694,
                             1.1936818095781789])
 
-        self.setSeed()
+        self.set_seed()
         actual = rayleigh(scale * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, rayleigh, bad_scale * 3)
@@ -1485,13 +1540,15 @@ def test_wald(self):
                             0.12450084820795027,
                             0.9096122728408238])
 
-        self.setSeed()
+        self.set_seed()
         actual = wald(mean * 3, scale)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, wald, bad_mean * 3, scale)
         assert_raises(ValueError, wald, mean * 3, bad_scale)
+        assert_raises(ValueError, random.wald, bad_mean * 3, scale)
+        assert_raises(ValueError, random.wald, mean * 3, bad_scale)
 
-        self.setSeed()
+        self.set_seed()
         actual = wald(mean, scale * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, wald, bad_mean, scale * 3)
@@ -1511,7 +1568,7 @@ def test_triangular(self):
                             2.0347400359389356,
                             2.0095991069536208])
 
-        self.setSeed()
+        self.set_seed()
         actual = triangular(left * 3, mode, right)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, triangular, bad_left_one * 3, mode, right)
@@ -1519,7 +1576,7 @@ def test_triangular(self):
         assert_raises(ValueError, triangular, bad_left_two * 3, bad_mode_two,
                       right)
 
-        self.setSeed()
+        self.set_seed()
         actual = triangular(left, mode * 3, right)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, triangular, bad_left_one, mode * 3, right)
@@ -1527,7 +1584,7 @@ def test_triangular(self):
         assert_raises(ValueError, triangular, bad_left_two, bad_mode_two * 3,
                       right)
 
-        self.setSeed()
+        self.set_seed()
         actual = triangular(left, mode, right * 3)
         assert_array_almost_equal(actual, desired, decimal=14)
         assert_raises(ValueError, triangular, bad_left_one, mode, right * 3)
@@ -1548,14 +1605,14 @@ def test_binomial(self):
         binom = random.binomial
         desired = np.array([1, 1, 1])
 
-        self.setSeed()
+        self.set_seed()
         actual = binom(n * 3, p)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, binom, bad_n * 3, p)
         assert_raises(ValueError, binom, n * 3, bad_p_one)
         assert_raises(ValueError, binom, n * 3, bad_p_two)
 
-        self.setSeed()
+        self.set_seed()
         actual = binom(n, p * 3)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, binom, bad_n, p * 3)
@@ -1571,14 +1628,14 @@ def test_negative_binomial(self):
         neg_binom = random.negative_binomial
         desired = np.array([1, 0, 1])
 
-        self.setSeed()
+        self.set_seed()
         actual = neg_binom(n * 3, p)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, neg_binom, bad_n * 3, p)
         assert_raises(ValueError, neg_binom, n * 3, bad_p_one)
         assert_raises(ValueError, neg_binom, n * 3, bad_p_two)
 
-        self.setSeed()
+        self.set_seed()
         actual = neg_binom(n, p * 3)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, neg_binom, bad_n, p * 3)
@@ -1594,7 +1651,7 @@ def test_poisson(self):
         poisson = random.poisson
         desired = np.array([1, 1, 0])
 
-        self.setSeed()
+        self.set_seed()
         actual = poisson(lam * 3)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, poisson, bad_lam_one * 3)
@@ -1606,7 +1663,7 @@ def test_zipf(self):
         zipf = random.zipf
         desired = np.array([2, 2, 1])
 
-        self.setSeed()
+        self.set_seed()
         actual = zipf(a * 3)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, zipf, bad_a * 3)
@@ -1621,7 +1678,7 @@ def test_geometric(self):
         geom = random.geometric
         desired = np.array([2, 2, 2])
 
-        self.setSeed()
+        self.set_seed()
         actual = geom(p * 3)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, geom, bad_p_one * 3)
@@ -1638,7 +1695,7 @@ def test_hypergeometric(self):
         hypergeom = random.hypergeometric
         desired = np.array([1, 1, 1])
 
-        self.setSeed()
+        self.set_seed()
         actual = hypergeom(ngood * 3, nbad, nsample)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, hypergeom, bad_ngood * 3, nbad, nsample)
@@ -1646,7 +1703,7 @@ def test_hypergeometric(self):
         assert_raises(ValueError, hypergeom, ngood * 3, nbad, bad_nsample_one)
         assert_raises(ValueError, hypergeom, ngood * 3, nbad, bad_nsample_two)
 
-        self.setSeed()
+        self.set_seed()
         actual = hypergeom(ngood, nbad * 3, nsample)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, hypergeom, bad_ngood, nbad * 3, nsample)
@@ -1654,7 +1711,7 @@ def test_hypergeometric(self):
         assert_raises(ValueError, hypergeom, ngood, nbad * 3, bad_nsample_one)
         assert_raises(ValueError, hypergeom, ngood, nbad * 3, bad_nsample_two)
 
-        self.setSeed()
+        self.set_seed()
         actual = hypergeom(ngood, nbad, nsample * 3)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, hypergeom, bad_ngood, nbad, nsample * 3)
@@ -1662,6 +1719,9 @@ def test_hypergeometric(self):
         assert_raises(ValueError, hypergeom, ngood, nbad, bad_nsample_one * 3)
         assert_raises(ValueError, hypergeom, ngood, nbad, bad_nsample_two * 3)
 
+        assert_raises(ValueError, hypergeom, -1, 10, 20)
+        assert_raises(ValueError, hypergeom, 10, -1, 20)
+        assert_raises(ValueError, hypergeom, 10, 10, 0)
         assert_raises(ValueError, hypergeom, 10, 10, 25)
 
     def test_logseries(self):
@@ -1671,7 +1731,7 @@ def test_logseries(self):
         logseries = random.logseries
         desired = np.array([1, 1, 1])
 
-        self.setSeed()
+        self.set_seed()
         actual = logseries(p * 3)
         assert_array_equal(actual, desired)
         assert_raises(ValueError, logseries, bad_p_one * 3)
@@ -1708,16 +1768,19 @@ def check_function(self, function, sz):
     def test_normal(self):
         def gen_random(state, out):
             out[...] = state.normal(size=10000)
+
         self.check_function(gen_random, sz=(10000,))
 
     def test_exp(self):
         def gen_random(state, out):
             out[...] = state.exponential(scale=np.ones((100, 1000)))
+
         self.check_function(gen_random, sz=(100, 1000))
 
     def test_multinomial(self):
         def gen_random(state, out):
-            out[...] = state.multinomial(10, [1/6.]*6, size=10000)
+            out[...] = state.multinomial(10, [1 / 6.] * 6, size=10000)
+
         self.check_function(gen_random, sz=(10000, 6))
 
 
@@ -1775,24 +1838,6 @@ def test_two_arg_funcs(self):
             out = func(self.argOne, argTwo[0])
             assert_equal(out.shape, self.tgtShape)
 
-# TODO: Uncomment once randint can broadcast arguments
-#    def test_randint(self):
-#        itype = [bool, np.int8, np.uint8, np.int16, np.uint16,
-#                 np.int32, np.uint32, np.int64, np.uint64]
-#        func = random.randint
-#        high = np.array([1])
-#        low = np.array([0])
-#
-#        for dt in itype:
-#            out = func(low, high, dtype=dt)
-#            self.assert_equal(out.shape, self.tgtShape)
-#
-#            out = func(low[0], high, dtype=dt)
-#            self.assert_equal(out.shape, self.tgtShape)
-#
-#            out = func(low, high[0], dtype=dt)
-#            self.assert_equal(out.shape, self.tgtShape)
-
     def test_three_arg_funcs(self):
         funcs = [random.noncentral_f, random.triangular,
                  random.hypergeometric]