fix!: Remove t_default and require t passed explicitly in SINDY.fit() and SINDY.score() (#628)

Previously, `t_default` was a convenience argument for people making demonstration notebooks.
It didn't do anything unique.  But by existing in __init__ and not fit(), it broke the distinction between
data-dependent and data-independent arguments.  But at the end of the day, it was just a duplicated
feature.

Author: yb6599

examples/1_feature_overview/example.ipynb

@@ -170,8 +170,8 @@ def ignore_specific_warnings():
 # Since we used a uniform timestep when defining `t_train` we can tell set a default timestep to be used whenever `t` isn't passed in.
 
 # In[10]:
-model = ps.SINDy(t_default=dt)
-model.fit(x_train)
+model = ps.SINDy()
+model.fit(x_train, t=dt)
 model.print()
 
 # %% [markdown]
@@ -255,7 +255,7 @@ def f(x):
 for i in range(1, n_steps):
     x_train_map[i] = f(x_train_map[i - 1]) + eps * np.random.randn()
 model = ps.SINDy(discrete_time=True)
-model.fit(x_train_map)
+model.fit(x_train_map, t=1)
 
 model.print()
 
@@ -1219,7 +1219,7 @@ def u_fun(t):
 model = ps.SINDy(
     feature_library=lib, optimizer=opt, feature_names=["x", "r"], discrete_time=True
 )
-model.fit(xs_train, u=rs_train)
+model.fit(xs_train, u=rs_train, t=1)
 model.print()
 
 # %% [markdown]
@@ -1407,7 +1407,7 @@ def u_fun(t):
 )
 optimizer = ps.STLSQ(threshold=8, alpha=1e-3, normalize_columns=True)
 model = ps.SINDy(feature_library=generalized_library, optimizer=optimizer)
-model.fit(data, x_dot=data_dot)
+model.fit(data, x_dot=data_dot, t=1)
 
 # Note scale of phi is large so some coefficients >> 1
 # --> would want to rescale phi with eps_0 for a harder problem

examples/1_feature_overview/example.py

@@ -36,9 +36,6 @@ class SINDyEstimator(SINDy):
         Names for the input features (e.g. ``['x', 'y', 'z']``). If None, will use
         ``['x0', 'x1', ...]``.
 
-    t_default : float, optional (default 1)
-        Default value for the time step.
-
     discrete_time : boolean, optional (default False)
         If True, dynamical system is treated as a map. Rather than predicting
         derivatives, the right hand side functions step the system forward by
@@ -65,15 +62,13 @@ def __init__(
         feature_library=None,
         differentiation_method=None,
         feature_names=None,
-        t_default=1,
         discrete_time=False,
     ):
         super(SINDyEstimator, self).__init__(
             optimizer=optimizer,
             feature_library=feature_library,
             differentiation_method=differentiation_method,
             feature_names=feature_names,
-            t_default=t_default,
             discrete_time=discrete_time,
         )
         self._model = None
@@ -102,7 +97,6 @@ def fit(self, x, **kwargs):
             feature_library=self.model.steps[0][1],
             optimizer=self.model.steps[-1][1],
             feature_names=self.feature_names,
-            t_default=self.t_default,
             discrete_time=self.discrete_time,
             n_control_features_=self.n_control_features_,
         )
@@ -166,9 +160,6 @@ class SINDyModel(SINDy):
         Names for the input features (e.g. ``['x', 'y', 'z']``). If None, will use
         ``['x0', 'x1', ...]``.
 
-    t_default : float, optional (default 1)
-        Default value for the time step.
-
     discrete_time : boolean, optional (default False)
         If True, dynamical system is treated as a map. Rather than predicting
         derivatives, the right hand side functions step the system forward by
@@ -193,15 +184,13 @@ def __init__(
         feature_library,
         optimizer,
         feature_names=None,
-        t_default=1,
         discrete_time=False,
         n_control_features_=0,
     ):
         super(SINDyModel, self).__init__(
             feature_library=feature_library,
             optimizer=optimizer,
             feature_names=feature_names,
-            t_default=t_default,
             discrete_time=discrete_time,
         )
         self.n_control_features_ = n_control_features_

pysindy/deeptime/deeptime.py

@@ -199,9 +199,6 @@ class SINDy(_BaseSINDy):
         Names for the input features (e.g. ``['x', 'y', 'z']``). If None, will use
         ``['x0', 'x1', ...]``.
 
-    t_default : float, optional (default 1)
-        Default value for the time step.
-
     discrete_time : boolean, optional (default False)
         If True, dynamical system is treated as a map. Rather than predicting
         derivatives, the right hand side functions step the system forward by
@@ -291,7 +288,6 @@ def __init__(
         feature_library: Optional[BaseFeatureLibrary] = None,
         differentiation_method: Optional[BaseDifferentiation] = None,
         feature_names: Optional[list[str]] = None,
-        t_default: float = 1,
         discrete_time: bool = False,
     ):
         if optimizer is None:
@@ -303,19 +299,13 @@ def __init__(
         if differentiation_method is None:
             differentiation_method = FiniteDifference(axis=-2)
         self.differentiation_method = differentiation_method
-        if not isinstance(t_default, float) and not isinstance(t_default, int):
-            raise ValueError("t_default must be a positive number")
-        elif t_default <= 0:
-            raise ValueError("t_default must be a positive number")
-        else:
-            self.t_default = t_default
         self.feature_names = feature_names
         self.discrete_time = discrete_time
 
     def fit(
         self,
         x,
-        t=None,
+        t,
         x_dot=None,
         u=None,
     ):
@@ -329,16 +319,14 @@ def fit(
             x should be a list containing data for each trajectory. Individual
             trajectories may contain different numbers of samples.
 
-        t: float, numpy array of shape (n_samples,), or list of numpy arrays, optional \
-                (default None)
+        t: float, numpy array of shape (n_samples,), or list of numpy arrays
             If t is a float, it specifies the timestep between each sample.
             If array-like, it specifies the time at which each sample was
             collected.
             In this case the values in t must be strictly increasing.
             In the case of multi-trajectory training data, t may also be a list
             of arrays containing the collection times for each individual
             trajectory.
-            If None, the default time step ``t_default`` will be used.
 
         x_dot: array-like or list of array-like, shape (n_samples, n_input_features), \
                 optional (default None)
@@ -362,9 +350,6 @@ def fit(
         self: a fitted :class:`SINDy` instance
         """
 
-        if t is None:
-            t = self.t_default
-
         if not _check_multiple_trajectories(x, x_dot, u):
             x, t, x_dot, u = _adapt_to_multiple_trajectories(x, t, x_dot, u)
         x, x_dot, u = _comprehend_and_validate_inputs(
@@ -482,7 +467,7 @@ def print(self, lhs=None, precision=3, **kwargs):
                 names = f"{lhs[i]}"
             print(f"{names} = {eqn}", **kwargs)
 
-    def score(self, x, t=None, x_dot=None, u=None, metric=r2_score, **metric_kws):
+    def score(self, x, t, x_dot=None, u=None, metric=r2_score, **metric_kws):
         """
         Returns a score for the time derivative prediction produced by the model.
 
@@ -491,12 +476,10 @@ def score(self, x, t=None, x_dot=None, u=None, metric=r2_score, **metric_kws):
         x: array-like or list of array-like, shape (n_samples, n_input_features)
             Samples from which to make predictions.
 
-        t: float, numpy array of shape (n_samples,), or list of numpy arrays, optional \
-                (default None)
+        t: float, numpy array of shape (n_samples,), or list of numpy arrays
             Time step between samples or array of collection times. Optional,
             used to compute the time derivatives of the samples if x_dot is not
             provided.
-            If None, the default time step ``t_default`` will be used.
 
         x_dot: array-like or list of array-like, shape (n_samples, n_input_features), \
                 optional (default None)
@@ -527,9 +510,6 @@ def score(self, x, t=None, x_dot=None, u=None, metric=r2_score, **metric_kws):
             Metric function value for the model prediction of x_dot.
         """
 
-        if t is None:
-            t = self.t_default
-
         if not _check_multiple_trajectories(x, x_dot, u):
             x, t, x_dot, u = _adapt_to_multiple_trajectories(x, t, x_dot, u)
         x, x_dot, u = _comprehend_and_validate_inputs(
@@ -597,7 +577,7 @@ def _process_trajectories(self, x, t, x_dot):
                 )
         return x, x_dot
 
-    def differentiate(self, x, t=None):
+    def differentiate(self, x, t):
         """
         Apply the model's differentiation method
         (:code:`self.differentiation_method`) to data.
@@ -607,10 +587,8 @@ def differentiate(self, x, t=None):
         x: array-like or list of array-like, shape (n_samples, n_input_features)
             Data to be differentiated.
 
-        t: int, numpy array of shape (n_samples,), or list of numpy arrays, optional \
-                (default None)
+        t: int, numpy array of shape (n_samples,), or list of numpy arrays
             Time step between samples or array of collection times.
-            If None, the default time step ``t_default`` will be used.
 
         Returns
         -------
@@ -622,8 +600,7 @@ def differentiate(self, x, t=None):
             "SINDy.differentiate is deprecated.  "
             "Call the differentiation_method parameter"
         )
-        if t is None:
-            t = self.t_default
+
         if self.discrete_time:
             raise RuntimeError("No differentiation implemented for discrete time model")
         if not _check_multiple_trajectories(x, None, None):

pysindy/pysindy.py

@@ -643,22 +643,22 @@ def test_parameterized_library(diffuse_multiple_trajectories):
 
 
 # Helper function for testing PDE libraries
-def pde_library_helper(library, u):
+def pde_library_helper(library, u, t):
     base_opt = STLSQ(normalize_columns=True, alpha=1e-10, threshold=0)
     model = SINDy(optimizer=base_opt, feature_library=library)
-    model.fit(u)
+    model.fit(u, t)
     assert np.any(base_opt.coef_ != 0.0)
 
 
 def test_1D_pdes(data_1d_random_pde):
-    _, spatial_grid, u, _ = data_1d_random_pde
+    t, spatial_grid, u, _ = data_1d_random_pde
     pde_lib = PDELibrary(
         function_library=PolynomialLibrary(degree=2, include_bias=False),
         derivative_order=4,
         spatial_grid=spatial_grid,
         include_bias=True,
     )
-    pde_library_helper(pde_lib, u)
+    pde_library_helper(pde_lib, u, t)
 
 
 def test_1D_weak_pdes():
@@ -676,7 +676,7 @@ def test_1D_weak_pdes():
         H_xt=2,
         include_bias=True,
     )
-    pde_library_helper(pde_lib, u)
+    pde_library_helper(pde_lib, u, t)
 
 
 def test_sindypi_library(data_lorenz):

test/test_feature_library.py

@@ -436,7 +436,7 @@ def test_sr3_quadratic_library(params):
     # Test SR3
     opt = SR3(**params)
     model = SINDy(optimizer=opt, feature_library=sindy_library)
-    model.fit(x)
+    model.fit(x, t=1)
     check_is_fitted(model)
 
 
@@ -472,7 +472,7 @@ def test_constrained_sr3_quadratic_library(params):
     # Test constrained SR3 without constraints
     opt = ConstrainedSR3(**params)
     model = SINDy(optimizer=opt, feature_library=sindy_library)
-    model.fit(x)
+    model.fit(x, t=1)
     check_is_fitted(model)
 
     # rerun with identity constraints
@@ -487,7 +487,7 @@ def test_constrained_sr3_quadratic_library(params):
         constraint_lhs=constraint_matrix, constraint_rhs=constraint_rhs, **params
     )
     model = SINDy(optimizer=opt, feature_library=sindy_library)
-    model.fit(x)
+    model.fit(x, t=1)
     check_is_fitted(model)
     assert np.allclose((model.coefficients().flatten())[:p], 0.0)
 
@@ -1081,7 +1081,7 @@ def test_ssr_criteria(data_lorenz):
     x, t = data_lorenz
     opt = SSR(normalize_columns=True, criteria="model_residual", kappa=1e-3)
     model = SINDy(optimizer=opt)
-    model.fit(x)
+    model.fit(x, t)
     assert np.shape(opt.coef_) == (3, 10)