Made the fixes from Brian in the pull request, i.e. fixed the noteboo…

…ks to not just ignore all warnings (only particular ones), fixed the SR3 examples in the Example 1 notebook (they were being used with the corrupted outlier data so were giving poor results), changed the random seed in Example 3 so that the logistic system is correctly identified, and added some warnings and syntax fixes to the TrappingSR3 optimizer.

examples/11_SSR_FROLS_examples.ipynb

@@ -2,6 +2,7 @@
  "cells": [
   {
    "cell_type": "markdown",
+   "id": "f2c6fdd0",
    "metadata": {},
    "source": [
     "### SSR and FROLs (the greedy algorithms!) examples\n",
@@ -13,6 +14,7 @@
   },
   {
    "cell_type": "markdown",
+   "id": "2d117798",
    "metadata": {},
    "source": [
     "Stepwise sparse regression (SSR) solves the problem by iteratively truncating the smallest coefficient during the optimization. There are many ways one can decide to truncate terms. We implement two popular ways: (1) truncating the smallest coefficient at each iteration; (2) chopping each coefficient, computing N - 1 models, and then choosing the model with the lowest residual error."
@@ -21,6 +23,7 @@
   {
    "cell_type": "code",
    "execution_count": 1,
+   "id": "af60b9cd",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -34,7 +37,8 @@
     "\n",
     "# Ignore odeint warnings when model is unstable\n",
     "import warnings\n",
-    "warnings.filterwarnings(\"ignore\")\n",
+    "from scipy.integrate.odepack import ODEintWarning\n",
+    "warnings.filterwarnings(\"ignore\", category=ODEintWarning)\n",
     "\n",
     "# Seed the random number generators for reproducibility\n",
     "np.random.seed(100)\n",
@@ -49,6 +53,7 @@
   {
    "cell_type": "code",
    "execution_count": 2,
+   "id": "c6a2adb0",
    "metadata": {},
    "outputs": [
     {
@@ -91,6 +96,7 @@
   },
   {
    "cell_type": "markdown",
+   "id": "7b2a932b",
    "metadata": {},
    "source": [
     "### Define some functions for plotting performance as the greedy algorithms progress"
@@ -99,6 +105,7 @@
   {
    "cell_type": "code",
    "execution_count": 3,
+   "id": "7b57f140",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -177,6 +184,7 @@
   },
   {
    "cell_type": "markdown",
+   "id": "f5dd9011",
    "metadata": {},
    "source": [
     "### Note that the usage is a bit different because we save all the sparse models and we choose our favorite one afterwards. Below we show that we can track the MSE between the predicted and true derivative on a testing trajectory as the algorithm iterates, and then choose the model with the minimum MSE. "
@@ -185,6 +193,7 @@
   {
    "cell_type": "code",
    "execution_count": 4,
+   "id": "73ac8fb3",
    "metadata": {},
    "outputs": [
     {
@@ -259,6 +268,7 @@
   {
    "cell_type": "code",
    "execution_count": 5,
+   "id": "0eda8aed",
    "metadata": {},
    "outputs": [
     {
@@ -297,6 +307,7 @@
   },
   {
    "cell_type": "markdown",
+   "id": "ed87b7b7",
    "metadata": {},
    "source": [
     "Note that some of the frames do not have any red lines... that means the model for this iteration resulted in an unstable model when it was integrated forward. This often happens with nonsparse models!\n",
@@ -307,6 +318,7 @@
   {
    "cell_type": "code",
    "execution_count": 6,
+   "id": "7ff558fb",
    "metadata": {},
    "outputs": [
     {
@@ -331,6 +343,7 @@
   {
    "cell_type": "code",
    "execution_count": 7,
+   "id": "e9ff19e5",
    "metadata": {},
    "outputs": [
     {
@@ -399,6 +412,7 @@
   {
    "cell_type": "code",
    "execution_count": 8,
+   "id": "8c394722",
    "metadata": {},
    "outputs": [
     {
@@ -436,6 +450,7 @@
   },
   {
    "cell_type": "markdown",
+   "id": "dd5126c4",
    "metadata": {},
    "source": [
     "# FROLS greedy algorithm\n",
@@ -446,6 +461,7 @@
   {
    "cell_type": "code",
    "execution_count": 9,
+   "id": "74f00770",
    "metadata": {
     "scrolled": true
    },
@@ -473,6 +489,7 @@
   {
    "cell_type": "code",
    "execution_count": 10,
+   "id": "ed98bf66",
    "metadata": {},
    "outputs": [
     {
@@ -541,6 +558,7 @@
   {
    "cell_type": "code",
    "execution_count": 11,
+   "id": "9bda98c0",
    "metadata": {},
    "outputs": [
     {
@@ -578,6 +596,7 @@
   },
   {
    "cell_type": "markdown",
+   "id": "1f0ef282",
    "metadata": {},
    "source": [
     "### Let's compare all three methods as the noise steadily increases, cross-validated over 10 noise instantiations"
@@ -586,6 +605,7 @@
   {
    "cell_type": "code",
    "execution_count": 12,
+   "id": "3970b7fa",
    "metadata": {},
    "outputs": [
     {
@@ -662,12 +682,12 @@
     "plt.xlabel('% training noise', fontsize=18)\n",
     "plt.ylabel('Minimum test MSE', fontsize=18)\n",
     "plt.title('Comparison of highest performance SSR and FROLs models,\\n '\n",
-    "          'averaged over 10 models trained on different noise', fontsize=16)\n",
-    "plt.savefig('SSR_FROLS_comparison.png')"
+    "          'averaged over 10 models trained on different noise', fontsize=16)"
    ]
   },
   {
    "cell_type": "markdown",
+   "id": "aa3336b4",
    "metadata": {},
    "source": [
     "### Summary\n",
@@ -691,7 +711,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.7.4"
+   "version": "3.7.7"
   },
   "toc": {
    "base_numbering": 1,

examples/13_ensembling.ipynb

@@ -51,7 +51,9 @@
     "\n",
     "# Ignore integration and solver convergence warnings\n",
     "import warnings\n",
-    "warnings.filterwarnings(\"ignore\")\n",
+    "from scipy.integrate.odepack import ODEintWarning\n",
+    "warnings.simplefilter(\"ignore\", category=UserWarning)\n",
+    "warnings.simplefilter(\"ignore\", category=ODEintWarning)\n",
     "\n",
     "import pysindy as ps\n",
     "\n",
@@ -484,7 +486,7 @@
     {
      "data": {
       "text/plain": [
-       "<matplotlib.legend.Legend at 0x7f5be57c3b50>"
+       "<matplotlib.legend.Legend at 0x7f9163f3db10>"
       ]
      },
      "execution_count": 12,

examples/5_differentiation.ipynb

@@ -38,12 +38,12 @@
     "\n",
     "# ignore user warnings\n",
     "import warnings\n",
-    "warnings.filterwarnings(\"ignore\")\n",
+    "warnings.filterwarnings(\"ignore\", category=UserWarning)\n",
     "\n",
     "integrator_keywords = {}\n",
-    "integrator_keywords['rtol'] = 1e-15\n",
+    "integrator_keywords['rtol'] = 1e-12\n",
     "integrator_keywords['method'] = 'LSODA'\n",
-    "integrator_keywords['atol'] = 1e-15"
+    "integrator_keywords['atol'] = 1e-12"
    ]
   },
   {
@@ -468,11 +468,11 @@
     "        feature_names=input_features\n",
     "    )\n",
     "    \n",
-    "    model.fit(x_train)\n",
+    "    model.fit(x_train, quiet=True)\n",
     "    equations_clean[name] = model.equations()\n",
     "    coefficients_clean[name] = model.coefficients()\n",
     "    \n",
-    "    model.fit(x_train_noisy)\n",
+    "    model.fit(x_train_noisy, quiet=True)\n",
     "    equations_noisy[name] = model.equations()\n",
     "    coefficients_noisy[name] = model.coefficients()"
    ]
@@ -661,18 +661,7 @@
      "start_time": "2020-10-22T00:48:21.700494Z"
     }
    },
-   "outputs": [
-    {
-     "name": "stderr",
-     "output_type": "stream",
-     "text": [
-      "/media/akaptano/HITSI-Work/pysindy/pysindy/optimizers/stlsq.py:176: UserWarning: Sparsity parameter is too big (0.5) and eliminated all coefficients\n",
-      "  \"coefficients\".format(self.threshold)\n",
-      "/media/akaptano/HITSI-Work/pysindy/pysindy/optimizers/stlsq.py:176: UserWarning: Sparsity parameter is too big (0.5) and eliminated all coefficients\n",
-      "  \"coefficients\".format(self.threshold)\n"
-     ]
-    }
-   ],
+   "outputs": [],
    "source": [
     "equations_clean = {}\n",
     "equations_noisy = {}\n",
@@ -690,11 +679,11 @@
     "        feature_names=input_features\n",
     "    )\n",
     "    \n",
-    "    model.fit(x_train)\n",
+    "    model.fit(x_train, quiet=True)\n",
     "    equations_clean[name] = model.equations()\n",
     "    coefficients_clean[name] = model.coefficients()\n",
     "    \n",
-    "    model.fit(x_train_noisy)\n",
+    "    model.fit(x_train_noisy, quiet=True)\n",
     "    equations_noisy[name] = model.equations()\n",
     "    coefficients_noisy[name] = model.coefficients()"
    ]

pysindy/__init__.py

@@ -4,7 +4,6 @@
 try:
     __version__ = get_distribution(__name__).version
 except DistributionNotFound:
-    __version__ = "v1.4.3"  # required for sphinx-build to compile locally
     pass
 
 from . import differentiation

pysindy/feature_library/base.py

@@ -355,7 +355,7 @@ def _name_combinations(self, lib_i, lib_j):
         lib_full = []
         for i in range(len(lib_i)):
             for j in range(len(lib_j)):
-                lib_full.append(lib_i[i] + lib_j[j])
+                lib_full.append(lib_i[i] + " " + lib_j[j])
         return lib_full
 
     def _set_inputs_per_library(self, inputs_per_library):

pysindy/optimizers/constrained_sr3.py

@@ -182,6 +182,9 @@ def __init__(
             self.unbias = False
             self.constraint_order = constraint_order
 
+        if inequality_constraints:
+            self.max_iter = max(10000, max_iter)  # max iterations for CVXPY
+
         if inequality_constraints and not self.use_constraints:
             raise ValueError(
                 "Use of inequality constraints requires constraint_lhs and "
@@ -240,7 +243,7 @@ def _update_coef_cvxpy(self, x, y, coef_sparse):
 
         # default solver is OSQP here but switches to ECOS for L2
         try:
-            prob.solve(max_iter=50000, eps_abs=self.tol, eps_rel=self.tol)
+            prob.solve(max_iter=self.max_iter, eps_abs=self.tol, eps_rel=self.tol)
         # Annoying error coming from L2 norm switching to use the ECOS
         # solver, which uses "max_iters" instead of "max_iter", and
         # similar semantic changes for the other variables.

pysindy/optimizers/trapping_sr3.py

@@ -338,7 +338,7 @@ def _set_Ptensors(self, r):
                             if (j == k) and (n == N - r + j) and (i == kk):
                                 PQ_tensor[i, j, k, kk, n] = 1.0
                             if (j != k) and (n == r + j + k - 1) and (i == kk):
-                                PQ_tensor[i, j, k, kk, n] = 1.0 / 2.0
+                                PQ_tensor[i, j, k, kk, n] = 1 / 2
 
         return PL_tensor_unsym, PL_tensor, PQ_tensor
 
@@ -357,6 +357,10 @@ def _check_P_matrix(self, r, n_features, N):
         # If these tensors are not passed, or incorrect shape, assume zeros
         if self.PQ_ is None:
             self.PQ_ = np.zeros((r, r, r, r, n_features))
+            warnings.warn(
+                "The PQ tensor (a requirement for the stability promotion) was"
+                " not set, so setting this tensor to all zeros. "
+            )
         elif (self.PQ_).shape != (r, r, r, r, n_features) and (self.PQ_).shape != (
             r,
             r,
@@ -365,15 +369,31 @@ def _check_P_matrix(self, r, n_features, N):
             N,
         ):
             self.PQ_ = np.zeros((r, r, r, r, n_features))
+            warnings.warn(
+                "The PQ tensor (a requirement for the stability promotion) was"
+                " initialized with incorrect dimensions, "
+                "so setting this tensor to all zeros "
+                "(with the correct dimensions). "
+            )
         if self.PL_ is None:
             self.PL_ = np.zeros((r, r, r, n_features))
+            warnings.warn(
+                "The PL tensor (a requirement for the stability promotion) was"
+                " not set, so setting this tensor to all zeros. "
+            )
         elif (self.PL_).shape != (r, r, r, n_features) and (self.PL_).shape != (
             r,
             r,
             r,
             N,
         ):
             self.PL_ = np.zeros((r, r, r, n_features))
+            warnings.warn(
+                "The PL tensor (a requirement for the stability promotion) was"
+                " initialized with incorrect dimensions, "
+                "so setting this tensor to all zeros "
+                "(with the correct dimensions). "
+            )
 
         # Check if the tensor symmetries are properly defined
         if self._bad_PL(self.PL_):
@@ -552,7 +572,7 @@ def _solve_direct_cvxpy(self, r, N, x_expanded, y, Pmatrix, coef_prev):
         cost = cp.sum_squares(x_expanded @ xi - y.flatten())
         if self.thresholder.lower() == "l1":
             cost = cost + self.threshold * cp.norm1(xi)
-        elif self.thresholder.lower() == "weighted_l1" and np.any(self.thresholds != 0):
+        elif self.thresholder.lower() == "weighted_l1":
             cost = cost + cp.norm1(np.ravel(self.thresholds) @ xi)
         elif self.thresholder.lower() == "l2":
             cost = cost + self.threshold * cp.norm2(xi) ** 2
@@ -590,7 +610,7 @@ def _solve_direct_cvxpy(self, r, N, x_expanded, y, Pmatrix, coef_prev):
             return None, None
         coef_sparse = (xi.value).reshape(coef_prev.shape)
 
-        if np.all(self.PL_ == 0.0) and np.all(self.PQ_ == 0.0):
+        if np.all(self.PL_) and np.all(self.PQ_):
             return np.zeros(r), coef_sparse  # no optimization over m
         else:
             m_cp = cp.Variable(r)