time delays not working correctly

[moonlians]

tutorial_koopman_hankel_dmdc_for_vdp_system.ipynb

%matplotlib inline
import pykoopman as pk
from pykoopman.common.examples import vdp_osc, rk4, square_wave  # required for example system
import matplotlib.pyplot as plt
import numpy as np
import numpy.random as rnd
np.random.seed(42)  # for reproducibility

import warnings
warnings.filterwarnings('ignore')

n_states = 2 # Number of states
n_inputs = 1 # Number of control inputs
dT = 0.01    # Timestep
n_traj = 200  # Number of trajectories
n_int = 1000  # Integration length
# Time vector
t = np.arange(0, n_int*dT, dT)

# Uniform random distributed forcing in [-1, 1]
u = 2*rnd.random([n_int, n_traj])-1

# Uniform distribution of initial conditions
x = 2*rnd.random([n_states, n_traj])-1

# Init
X = np.zeros((n_states, n_int*n_traj))
Y = np.zeros((n_states, n_int*n_traj))
U = np.zeros((n_inputs, n_int*n_traj))

# Integrate
for step in range(n_int):
    y = rk4(0, x, u[step, :], dT, vdp_osc)
    X[:, (step)*n_traj:(step+1)*n_traj] = x
    Y[:, (step)*n_traj:(step+1)*n_traj] = y
    U[:, (step)*n_traj:(step+1)*n_traj] = u[step, :]
    x = y
# Visualize first 100 steps of the training data
fig, axs = plt.subplots(1, 1, tight_layout=True, figsize=(12, 4))
for traj_idx in range(n_traj):
    x = X[:, traj_idx::n_traj]
    axs.plot(t[0:100], x[1, 0:100], 'k')
axs.set(
        ylabel=r'$x_2$',
        xlabel=r'$t$')

EDMDc = pk.regression.DMDc()

n_delays = 9
obs = pk.observables.TimeDelay(n_delays=n_delays)

model = pk.Koopman(observables=obs, regressor=EDMDc)
model.fit(X.T, y=Y.T, u=U[:,n_delays:].T)
Traceback (most recent call last):

  Cell In[291], line 53
    model.fit(X.T, y=Y.T, u=U[:,n_delays:].T)

  File ~\anaconda3\envs\DMD\lib\site-packages\pykoopman\koopman.py:153 in fit
    self.model.fit(x, y, regressor__u=u, regressor__dt=dt)

  File ~\anaconda3\envs\DMD\lib\site-packages\sklearn\pipeline.py:394 in fit
    self._final_estimator.fit(Xt, y, **fit_params_last_step)

  File ~\anaconda3\envs\DMD\lib\site-packages\sklearn\compose\_target.py:229 in fit
    self._fit_transformer(y_2d)

  File ~\anaconda3\envs\DMD\lib\site-packages\sklearn\compose\_target.py:175 in _fit_transformer
    self.transformer_.fit(y)

  File ~\anaconda3\envs\DMD\lib\site-packages\sklearn\preprocessing\_function_transformer.py:165 in fit
    self._check_inverse_transform(X)

  File ~\anaconda3\envs\DMD\lib\site-packages\sklearn\preprocessing\_function_transformer.py:136 in _check_inverse_transform
    if not _allclose_dense_sparse(X[idx_selected], X_round_trip):

  File ~\anaconda3\envs\DMD\lib\site-packages\sklearn\utils\validation.py:1607 in _allclose_dense_sparse
    return np.allclose(x, y, rtol=rtol, atol=atol)

  File <__array_function__ internals>:180 in allclose

  File ~\anaconda3\envs\DMD\lib\site-packages\numpy\core\numeric.py:2265 in allclose
    res = all(isclose(a, b, rtol=rtol, atol=atol, equal_nan=equal_nan))

  File <__array_function__ internals>:180 in isclose

  File ~\anaconda3\envs\DMD\lib\site-packages\numpy\core\numeric.py:2375 in isclose
    return within_tol(x, y, atol, rtol)

  File ~\anaconda3\envs\DMD\lib\site-packages\numpy\core\numeric.py:2356 in within_tol
    return less_equal(abs(x-y), atol + rtol * abs(y))

ValueError: operands could not be broadcast together with shapes (100,2) (91,2) 

I tried also normal EDMD with time delay as an observable and it errored with the same message.

EDMD = pk.regression.EDMD()
model = pk.Koopman(observables=obs, regressor=EDMD)
model.fit(X.T, y=Y.T)

[pswpswpsw]

It is because the pykoopman you imported is the one from PyPi, which is a very old version. I just uploaded the 1.0.1 version. You should be able to find the example working now.

[moonlians]

Thanks! I uninstalled and reinstalled, tried both from git clone and pip install pykoopman.
But it for some reason is not loading correctly:
Traceback (most recent call last):

Cell In[891], line 1
EDMDc = pk.regression.DMDc()

AttributeError: module 'pykoopman' has no attribute 'regression'

I definitely have "regression" folder inside the git cloned folder...

EDIT: it works as from pykoopman.regression import DMDc, so I can do that.

[moonlians]

There seems to also be a problem instituting the SVD rank:

from pykoopman.regression import EDMDc
EDMDc = DMDc(svd_rank=20)

n_delays = 71
obs = TimeDelay(n_delays=n_delays)

model = pk.Koopman(observables=obs, regressor=EDMDc)
model.fit(X1.T, y=X2.T, u=U[:,n_delays:len(X)].T) 

gives the error:

model.fit(X1.T, y=X2.T, u=U[:,n_delays:len(X)].T) 
Traceback (most recent call last):

  Cell In[3541], line 1
    model.fit(X1.T, y=X2.T, u=U[:,n_delays:len(X)].T)

  File ~\anaconda3\envs\DMD\lib\site-packages\pykoopman\koopman.py:170 in fit
    self._pipeline.fit(x, y, regressor__u=u, regressor__dt=dt)

  File ~\anaconda3\envs\DMD\lib\site-packages\sklearn\pipeline.py:394 in fit
    self._final_estimator.fit(Xt, y, **fit_params_last_step)

  File ~\anaconda3\envs\DMD\lib\site-packages\pykoopman\regression\_base_ensemble.py:107 in fit
    self.regressor_.fit(X, y_trans, **fit_params)

  File ~\anaconda3\envs\DMD\lib\site-packages\pykoopman\regression\_dmdc.py:146 in fit
    self._fit_unknown_B(X1, X2, C, r, rout)

  File ~\anaconda3\envs\DMD\lib\site-packages\pykoopman\regression\_dmdc.py:175 in _fit_unknown_B
    assert rout <= r

AssertionError

It works fine when I don't include the svd_rank parameter.

[pswpswpsw]

can you have a complete self-contained code here? including the data

[pswpswpsw]

please take a look at the cell [5] in https://github.com/dynamicslab/pykoopman/blob/master/docs/tutorial_koopman_hankel_dmdc_for_vdp_system.ipynb

So basically, if you are implementing EDMDc with time delay, set svd_output_rank or svd_rank less than the full rank might lead to worse performance.

In fact, this is not a special case. Here is an example in my JFM paper:

• Looking at figure 12 in the following paper, when using EDMD with monomials, choosing a svd rank in the least square does not help.

Reference

• Sparsity-promoting algorithms for the discovery of informative Koopman-invariant subspaces
  S Pan, N Arnold-Medabalimi, K Duraisamy
  Journal of Fluid Mechanics 917, A18

[moonlians]

Thanks, that's really helpful!

[moonlians]

So basically, if you are implementing EDMDc with time delay, set svd_output_rank or svd_rank less than the full rank might lead to worse performance.

It seems like when there is noise in the control signal, it unfortunately leads to overfitting