pysindy.optimizers.TrappingSR3

class pysindy.optimizers.TrappingSR3(*, _n_tgts: int = None, _include_bias: bool = False, _interaction_only: bool = False, method: str = 'global', eta: float | None = None, eps_solver: float = 1e-07, alpha: float | None = None, beta: float | None = None, mod_matrix: ndarray[tuple[Any, ...], dtype[_ScalarT]] | None = None, alpha_A: float | None = None, alpha_m: float | None = None, gamma: float = -0.1, tol_m: float = 1e-05, regularizer: str = 'l1', m0: ndarray[tuple[Any, ...], dtype[_ScalarT]] | None = None, A0: ndarray[tuple[Any, ...], dtype[_ScalarT]] | None = None, **kwargs)[source]

Generalized trapping variant of sparse relaxed regularized regression. This optimizer can be used to identify quadratically nonlinear systems with either a-priori globally or locally stable (bounded) solutions.

This optimizer can be used to minimize five different objective functions:

\[0.5\|y-Xw\|^2_2 + \lambda \times R(w) + 0.5\|Pw-A\|^2_2/\eta + \delta_0(Cw-d) + \delta_{\Lambda}(A) + \alpha \|Qijk\| + \beta \|Q_{ijk} + Q_{jik} + Q_{kij}\|\]

where $R(w)$ is a regularization function, C is a constraint matrix detailing affine constraints on the model coefficients, A is a proxy for the quadratic contributions to the energy evolution, and $Q_{ijk}$ are the quadratic coefficients in the model. For provably globally bounded solutions, use $\alpha >> 1$, $\beta >> 1$, and equality constraints. For maximizing the local stability radius of the model one has the choice to do this by (1) minimizing the values in $Q_{ijk}$, (2) promoting models with skew-symmetrix $Q_{ijk}$ coefficients, or (3) using inequality constraints for skew-symmetry in $Q_{ijk}$.

See the following references for more details:

Kaptanoglu, Alan A., et al. “Promoting global stability in data-driven models of quadratic nonlinear dynamics.” arXiv preprint arXiv:2105.01843 (2021).

Parameters:

eta – Determines the strength of the stability term $||Pw-A||^2$ in the optimization. The default value is very large so that the algorithm default is to ignore the stability term. In this limit, this should be approximately equivalent to the ConstrainedSR3 method.
eps_solver – If reg_weight_lam != 0, this specifies the error tolerance in the CVXPY (OSQP) solve. Default 1.0e-7
alpha – Determines the strength of the local stability term $||Qijk||^2$ in the optimization. The default value (1e20) is very large so that the algorithm default is to ignore this term.
beta – Determines the strength of the local stability term $||Qijk + Qjik + Qkij||^2$ in the optimization. The default value is very large so that the algorithm default is to ignore this term.
mod_matrix – Lyapunov matrix. Trapping theorems apply to energy $\propto \dot y \cdot y$, but also to any $\propto \dot y P \cdot y$ for Lyapunov matrix $P$. Defaults to the identity matrix.
alpha_A – Determines the step size in the prox-gradient descent over A. For convergence, need alpha_A <= eta, so default alpha_A = eta is used.
alpha_m – Determines the step size in the prox-gradient descent over m. For convergence, need alpha_m <= eta / ||w^T * PQ^T * PQ * w||. Typically 0.01 * eta <= alpha_m <= 0.1 * eta. (default eta * 0.1)
gamma – Determines the negative interval that matrix A is projected onto. For most applications gamma = 0.1 - 1.0 works pretty well.
tol_m – Tolerance used for determining convergence of the optimization algorithm over m.
regularizer – Regularization function to use. For current trapping SINDy, only the L1 and L2 norms are implemented. Note that other convex norms could be straightforwardly implemented, but L0 requires reformulation because of nonconvexity. (default ‘L1’)
m0 – Initial guess for trap center in the optimization. Default None initializes vector elements randomly in [-1, 1]. shape (n_targets)
A0 – Initial guess for vector A in the optimization. Shape (n_targets, n_targets) Default None, meaning A is initialized as A = diag(gamma).

Attributes:

A_history_ (list) – History of the auxiliary variable A that approximates diag(PW).
m_history_ (list) – History of the shift vector m that determines the origin of the trapping region.
PW_history_ (list) – History of PW = A^S, the quantity we are attempting to make negative definite.
PWeigs_history_ (list) – History of diag(PW), a list of the eigenvalues of A^S at each iteration. Tracking this allows us to ascertain if A^S is indeed being pulled towards the space of negative definite matrices.
PL_unsym_ (np.ndarray, shape (n_targets, n_targets, n_targets, n_features)) – Unsymmetrized linear coefficient part of the P matrix in :math:||Pw - A||^2
PL_ (np.ndarray, shape (n_targets, n_targets, n_targets, n_features)) – Linear coefficient part of the P matrix in :math:||Pw - A||^2
PQ_ (np.ndarray, shape (n_targets, n_targets,) – n_targets, n_targets, n_features) Quadratic coefficient part of the P matrix in :math:||Pw - A||^2
PT_ (np.ndarray, shape (n_targets, n_targets,) – n_targets, n_targets, n_features) Transpose of 1st dimension and 2nd dimension of quadratic coefficient part of the P matrix in :math:||Pw - A||^2
objective_history_ (list) – History of the value of the objective at each step. Note that the trapping SINDy problem is nonconvex, meaning that this value may increase and decrease as the algorithm works.

Examples

>>> import numpy as np
>>> from scipy.integrate import odeint
>>> from pysindy import SINDy
>>> from pysindy.optimizers import TrappingSR3
>>> lorenz = lambda z,t : [10*(z[1] - z[0]),
>>>                        z[0]*(28 - z[2]) - z[1],
>>>                        z[0]*z[1] - 8/3*z[2]]
>>> t = np.arange(0,2,.002)
>>> x = odeint(lorenz, [-8,8,27], t)
>>> opt = TrappingSR3(reg_weight_lam=0.1)
>>> model = SINDy(optimizer=opt)
>>> model.fit(x, t=t[1]-t[0])
>>> model.print()
x0' = -10.004 1 + 10.004 x0
x1' = 27.994 1 + -0.993 x0 + -1.000 1 x1
x2' = -2.662 x1 + 1.000 1 x0

Methods

`set_fit_request`	Configure whether metadata should be requested to be passed to the `fit` method.
`set_params`	Set the parameters of this estimator.
`set_score_request`	Configure whether metadata should be requested to be passed to the `score` method.

Attributes

`max_iter`
`normalize_columns`
`initial_guess`
`copy_X`
`unbias`
`coef_`
`intercept_`

set_fit_request(*, sample_weight: bool | None | str = '$UNCHANGED$', x_: bool | None | str = '$UNCHANGED$') → TrappingSR3

Configure whether metadata should be requested to be passed to the fit method.

Note that this method is only relevant when this estimator is used as a sub-estimator within a meta-estimator and metadata routing is enabled with enable_metadata_routing=True (see sklearn.set_config). Please check the User Guide on how the routing mechanism works.

The options for each parameter are:

True: metadata is requested, and passed to fit if provided. The request is ignored if metadata is not provided.
False: metadata is not requested and the meta-estimator will not pass it to fit.
None: metadata is not requested, and the meta-estimator will raise an error if the user provides it.
str: metadata should be passed to the meta-estimator with this given alias instead of the original name.

The default (sklearn.utils.metadata_routing.UNCHANGED) retains the existing request. This allows you to change the request for some parameters and not others.

Added in version 1.3.

Parameters:

sample_weight (str, True, False, or None, default=sklearn.utils.metadata_routing.UNCHANGED) – Metadata routing for sample_weight parameter in fit.
x (str, True, False, or None, default=sklearn.utils.metadata_routing.UNCHANGED) – Metadata routing for x_ parameter in fit.

Returns:

self – The updated object.

Return type:

object

set_params(**kwargs)[source]

Set the parameters of this estimator.

The method works on simple estimators as well as on nested objects (such as Pipeline). The latter have parameters of the form <component>__<parameter> so that it’s possible to update each component of a nested object.

Parameters:: **params (dict) – Estimator parameters.
Returns:: self – Estimator instance.
Return type:: estimator instance

set_score_request(*, sample_weight: bool | None | str = '$UNCHANGED$') → TrappingSR3

Configure whether metadata should be requested to be passed to the score method.

Note that this method is only relevant when this estimator is used as a sub-estimator within a meta-estimator and metadata routing is enabled with enable_metadata_routing=True (see sklearn.set_config). Please check the User Guide on how the routing mechanism works.

The options for each parameter are:

True: metadata is requested, and passed to score if provided. The request is ignored if metadata is not provided.
False: metadata is not requested and the meta-estimator will not pass it to score.
None: metadata is not requested, and the meta-estimator will raise an error if the user provides it.
str: metadata should be passed to the meta-estimator with this given alias instead of the original name.

The default (sklearn.utils.metadata_routing.UNCHANGED) retains the existing request. This allows you to change the request for some parameters and not others.

Added in version 1.3.

Parameters:: sample_weight (str, True, False, or None, default=sklearn.utils.metadata_routing.UNCHANGED) – Metadata routing for sample_weight parameter in score.
Returns:: self – The updated object.
Return type:: object