[numpy.md] Update np.random → Generator API (#549)

Author: Chihiro2000GitHub

lectures/numpy.md

@@ -864,7 +864,8 @@ a
 Mutability leads to the following behavior (which can be shocking to MATLAB programmers...)
 
 ```{code-cell} python3
-a = np.random.randn(3)
+rng = np.random.default_rng()
+a = rng.standard_normal(3)
 a
 ```
 
@@ -894,7 +895,7 @@ It is of course possible to make `b` an independent copy of `a` when required.
 This can be done using `np.copy`
 
 ```{code-cell} python3
-a = np.random.randn(3)
+a = rng.standard_normal(3)
 a
 ```
 
@@ -972,7 +973,7 @@ def f(x):
 The NumPy function `np.where` provides a vectorized alternative:
 
 ```{code-cell} python3
-x = np.random.randn(4)
+x = rng.standard_normal(4)
 x
 ```
 
@@ -1049,11 +1050,12 @@ z[z > 3]
 NumPy provides some additional functionality related to scientific programming
 through its sub-packages.
 
-We've already seen how we can generate random variables using np.random
+We've already seen how we can generate random variables using NumPy's 
+[random `Generator`](https://numpy.org/doc/stable/reference/random/generator.html#random-generator).
 
 ```{code-cell} python3
-z = np.random.randn(10000)  # Generate standard normals
-y = np.random.binomial(10, 0.5, size=1000)    # 1,000 draws from Bin(10, 0.5)
+z = rng.standard_normal(10000)  # Generate standard normals
+y = rng.binomial(10, 0.5, size=1000)    # 1,000 draws from Bin(10, 0.5)
 y.mean()
 ```
 
@@ -1099,7 +1101,7 @@ It takes a few seconds to run.
 n = 20
 m = 1000
 for i in range(n):
-    X = np.random.randn(m, m)
+    X = rng.standard_normal((m, m))
     λ = np.linalg.eigvals(X)
 ```
 
@@ -1240,28 +1242,28 @@ Here's our first pass at a solution:
 
 ```{code-cell} python3
 from numpy import cumsum
-from numpy.random import uniform
 
 class DiscreteRV:
     """
     Generates an array of draws from a discrete random variable with vector of
     probabilities given by q.
     """
 
-    def __init__(self, q):
+    def __init__(self, q, seed=None):
         """
         The argument q is a NumPy array, or array like, nonnegative and sums
         to 1
         """
         self.q = q
         self.Q = cumsum(q)
+        self.rng = np.random.default_rng(seed)
 
     def draw(self, k=1):
         """
         Returns k draws from q. For each such draw, the value i is returned
         with probability q[i].
         """
-        return self.Q.searchsorted(uniform(0, 1, size=k))
+        return self.Q.searchsorted(self.rng.uniform(0, 1, size=k))
 ```
 
 The logic is not obvious, but if you take your time and read it slowly,
@@ -1390,7 +1392,8 @@ Here's an example of usage
 
 ```{code-cell} python3
 fig, ax = plt.subplots()
-X = np.random.randn(1000)
+rng = np.random.default_rng()
+X = rng.standard_normal(1000)
 F = ECDF(X)
 F.plot(ax)
 ```
@@ -1411,9 +1414,9 @@ In this exercise, try to use `for` loops to replicate the result of the followin
 
 ```{code-cell} python3
 
-np.random.seed(123)
-x = np.random.randn(4, 4)
-y = np.random.randn(4)
+rng = np.random.default_rng(123)
+x = rng.standard_normal((4, 4))
+y = rng.standard_normal(4)
 A = x / y
 ```
 
@@ -1441,9 +1444,9 @@ Now we can import the quantecon package.
 
 ```{code-cell} python3
 
-np.random.seed(123)
-x = np.random.randn(1000, 100, 100)
-y = np.random.randn(100)
+rng = np.random.default_rng(123)
+x = rng.standard_normal((1000, 100, 100))
+y = rng.standard_normal(100)
 
 with qe.Timer("Broadcasting operation"):
     B = x / y
@@ -1469,9 +1472,9 @@ print(B)
 **Part 1 Solution**
 
 ```{code-cell} python3
-np.random.seed(123)
-x = np.random.randn(4, 4)
-y = np.random.randn(4)
+rng = np.random.default_rng(123)
+x = rng.standard_normal((4, 4))
+y = rng.standard_normal(4)
 
 C = np.empty_like(x)
 n = len(x)
@@ -1500,9 +1503,9 @@ print(np.array_equal(A, C))
 
 ```{code-cell} python3
 
-np.random.seed(123)
-x = np.random.randn(1000, 100, 100)
-y = np.random.randn(100)
+rng = np.random.default_rng(123)
+x = rng.standard_normal((1000, 100, 100))
+y = rng.standard_normal(100)
 
 with qe.Timer("For loop operation"):
     D = np.empty_like(x)