diff --git a/pyproject.toml b/pyproject.toml
index c1ab97e23f9366426b5fc9b4e52e4c06990fed52..3cfd42e4891039f95f9c428606dcc268db3ffda2 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -100,7 +100,7 @@ precision = 2
 exclude = ["examples"]
 
 [tool.ruff.lint]
-select = ["E4", "E7", "E9", "F", "SIM", "B"]
+select = ["E4", "E7", "E9", "F", "SIM", "B", "NPY"]
 ignore = ["F403", "B008", "B021", "B006"]
 
 [tool.typos]
diff --git a/test/unit/test_list_schedulers.py b/test/unit/test_list_schedulers.py
index 52d29a130c9dad0dff2c6c3ec2ea6784502640ed..b67147d9e5d297b1140eac72fa790d902f50eaec 100644
--- a/test/unit/test_list_schedulers.py
+++ b/test/unit/test_list_schedulers.py
@@ -882,7 +882,7 @@ class TestHybridScheduler:
         assert schedule.schedule_time == 16
 
         # validate regenerated sfg with random 2x2 real s.p.d. matrix
-        A = np.random.rand(2, 2)
+        A = np.random.default_rng().random((2, 2))
         A = np.dot(A, A.T)
         A_inv = np.linalg.inv(A)
         input_signals = []
@@ -1985,7 +1985,7 @@ def _validate_recreated_sfg_ldlt_matrix_inverse(
         delays = [0 for i in range(num_of_outputs)]
 
     # random real s.p.d matrix
-    A = np.random.rand(N, N)
+    A = np.random.default_rng().random((N, N))
     A = np.dot(A, A.T)
 
     # iterate through the upper diagonal and construct the input to the SFG
diff --git a/test/unit/test_sfg_generators.py b/test/unit/test_sfg_generators.py
index a51224bd121887976bc3879b4316c757ee999036..f427259e9cc6917a01f6e1a616f09dcc903b78dd 100644
--- a/test/unit/test_sfg_generators.py
+++ b/test/unit/test_sfg_generators.py
@@ -310,7 +310,7 @@ class TestDirectFormIIRType1:
 
         b, a = signal.butter(N, Wc, btype="lowpass", output="ba")
 
-        input_signal = np.random.randn(100)
+        input_signal = np.random.default_rng().standard_normal(100)
         reference_filter_output = signal.lfilter(b, a, input_signal)
 
         sfg = direct_form_1_iir(b, a, name="test iir direct form 1")
@@ -326,7 +326,7 @@ class TestDirectFormIIRType1:
 
         b, a = signal.butter(N, Wc, btype="lowpass", output="ba")
 
-        input_signal = np.random.randn(100)
+        input_signal = np.random.default_rng().standard_normal(100)
         reference_filter_output = signal.lfilter(b, a, input_signal)
 
         sfg = direct_form_1_iir(b, a, name="test iir direct form 1")
@@ -343,7 +343,7 @@ class TestDirectFormIIRType1:
         b, a = signal.ellip(N, 0.1, 60, Wc, btype="low", analog=False)
         b, a = signal.butter(N, Wc, btype="lowpass", output="ba")
 
-        input_signal = np.random.randn(100)
+        input_signal = np.random.default_rng().standard_normal(100)
         reference_filter_output = signal.lfilter(b, a, input_signal)
 
         sfg = direct_form_1_iir(b, a, name="test iir direct form 1")
@@ -430,7 +430,7 @@ class TestDirectFormIIRType2:
 
         b, a = signal.butter(N, Wc, btype="lowpass", output="ba")
 
-        input_signal = np.random.randn(100)
+        input_signal = np.random.default_rng().standard_normal(100)
         reference_filter_output = signal.lfilter(b, a, input_signal)
 
         sfg = direct_form_2_iir(b, a, name="test iir direct form 1")
@@ -446,7 +446,7 @@ class TestDirectFormIIRType2:
 
         b, a = signal.butter(N, Wc, btype="lowpass", output="ba")
 
-        input_signal = np.random.randn(100)
+        input_signal = np.random.default_rng().standard_normal(100)
         reference_filter_output = signal.lfilter(b, a, input_signal)
 
         sfg = direct_form_2_iir(b, a, name="test iir direct form 1")
@@ -463,7 +463,7 @@ class TestDirectFormIIRType2:
         b, a = signal.ellip(N, 0.1, 60, Wc, btype="low", analog=False)
         b, a = signal.butter(N, Wc, btype="lowpass", output="ba")
 
-        input_signal = np.random.randn(100)
+        input_signal = np.random.default_rng().standard_normal(100)
         reference_filter_output = signal.lfilter(b, a, input_signal)
 
         sfg = direct_form_2_iir(b, a, name="test iir direct form 1")
@@ -804,7 +804,7 @@ class TestLdltMatrixInverse:
     #         assert np.isclose(actual_values[i], expected_values[i])
 
     def _generate_random_spd_matrix(self, N: int) -> np.ndarray:
-        A = np.random.rand(N, N)
+        A = np.random.default_rng().random((N, N))
         A = (A + A.T) / 2  # ensure symmetric
         min_eig = np.min(np.linalg.eigvals(A))
         A += (np.abs(min_eig) + 0.1) * np.eye(N)  # ensure positive definiteness