diff --git a/control/rlocus.py b/control/rlocus.py
index 631185cc6..281fed082 100644
--- a/control/rlocus.py
+++ b/control/rlocus.py
@@ -21,7 +21,7 @@
 import matplotlib.pyplot as plt
 import numpy as np
 import scipy.signal  # signal processing toolbox
-from numpy import array, imag, poly1d, real, row_stack, zeros_like
+from numpy import array, imag, poly1d, real, vstack, zeros_like
 
 from . import config
 from .exception import ControlMIMONotImplemented
@@ -421,7 +421,7 @@ def _RLFindRoots(nump, denp, kvect):
         curroots.sort()
         roots.append(curroots)
 
-    return row_stack(roots)
+    return vstack(roots)
 
 
 def _RLSortRoots(roots):
diff --git a/control/statesp.py b/control/statesp.py
index e14a8358a..0c2856b15 100644
--- a/control/statesp.py
+++ b/control/statesp.py
@@ -48,26 +48,27 @@
 """
 
 import math
+from copy import deepcopy
+from warnings import warn
+
 import numpy as np
-from numpy import any, asarray, concatenate, cos, delete, \
-    empty, exp, eye, isinf, ones, pad, sin, zeros, squeeze
-from numpy.random import rand, randn
-from numpy.linalg import solve, eigvals, matrix_rank
-from numpy.linalg.linalg import LinAlgError
 import scipy as sp
 import scipy.linalg
-from scipy.signal import cont2discrete
+from numpy import (any, asarray, concatenate, cos, delete, empty, exp, eye,
+                   isinf, ones, pad, sin, squeeze, zeros)
+from numpy.linalg import LinAlgError, eigvals, matrix_rank, solve
+from numpy.random import rand, randn
 from scipy.signal import StateSpace as signalStateSpace
-from warnings import warn
+from scipy.signal import cont2discrete
 
-from .exception import ControlSlycot, slycot_check, ControlMIMONotImplemented
+from . import config
+from .exception import ControlMIMONotImplemented, ControlSlycot, slycot_check
 from .frdata import FrequencyResponseData
+from .iosys import (InputOutputSystem, _process_dt_keyword,
+                    _process_iosys_keywords, _process_signal_list,
+                    common_timebase, isdtime, issiso)
 from .lti import LTI, _process_frequency_response
-from .iosys import InputOutputSystem, common_timebase, isdtime, issiso, \
-    _process_iosys_keywords, _process_dt_keyword, _process_signal_list
-from .nlsys import NonlinearIOSystem, InterconnectedSystem
-from . import config
-from copy import deepcopy
+from .nlsys import InterconnectedSystem, NonlinearIOSystem
 
 try:
     from slycot import ab13dd
@@ -2221,9 +2222,10 @@ def _convert_to_statespace(sys, use_prefix_suffix=False, method=None):
     by the calling function.
 
     """
-    from .xferfcn import TransferFunction
     import itertools
 
+    from .xferfcn import TransferFunction
+
     if isinstance(sys, StateSpace):
         return sys
 
diff --git a/control/tests/timeresp_test.py b/control/tests/timeresp_test.py
index fb21180b3..bdbbb3e89 100644
--- a/control/tests/timeresp_test.py
+++ b/control/tests/timeresp_test.py
@@ -173,15 +173,15 @@ def tsystem(self, request):
         # System Type 1 - Step response not stationary:  G(s)=1/s(s+1)
         siso_tf_type1 = TSys(TransferFunction(1, [1, 1, 0]))
         siso_tf_type1.step_info = {
-             'RiseTime': np.NaN,
-             'SettlingTime': np.NaN,
-             'SettlingMin': np.NaN,
-             'SettlingMax': np.NaN,
-             'Overshoot': np.NaN,
-             'Undershoot': np.NaN,
-             'Peak': np.Inf,
-             'PeakTime': np.Inf,
-             'SteadyStateValue': np.NaN}
+             'RiseTime': np.nan,
+             'SettlingTime': np.nan,
+             'SettlingMin': np.nan,
+             'SettlingMax': np.nan,
+             'Overshoot': np.nan,
+             'Undershoot': np.nan,
+             'Peak': np.inf,
+             'PeakTime': np.inf,
+             'SteadyStateValue': np.nan}
 
         # SISO under shoot response and positive final value
         # G(s)=(-s+1)/(s²+s+1)
diff --git a/control/timeresp.py b/control/timeresp.py
index 58207e88e..428baf230 100644
--- a/control/timeresp.py
+++ b/control/timeresp.py
@@ -1590,15 +1590,15 @@ def step_info(sysdata, T=None, T_num=None, yfinal=None, params=None,
             InfValue = InfValues[i, j]
             sgnInf = np.sign(InfValue.real)
 
-            rise_time: float = np.NaN
-            settling_time: float = np.NaN
-            settling_min: float = np.NaN
-            settling_max: float = np.NaN
-            peak_value: float = np.Inf
-            peak_time: float = np.Inf
-            undershoot: float = np.NaN
-            overshoot: float = np.NaN
-            steady_state_value: complex = np.NaN
+            rise_time: float = np.nan
+            settling_time: float = np.nan
+            settling_min: float = np.nan
+            settling_max: float = np.nan
+            peak_value: float = np.inf
+            peak_time: float = np.inf
+            undershoot: float = np.nan
+            overshoot: float = np.nan
+            steady_state_value: complex = np.nan
 
             if not np.isnan(InfValue) and not np.isinf(InfValue):
                 # RiseTime