A simulation of the discovered model diverges rapidly after first timestep and cannot capture the entire dynamics -- Super Unstable Model #330
Replies: 3 comments 2 replies
-
|
Thanks for the question! It looks like your data is either changing very quickly or your measurement noise is pretty high. If the latter, you could try smoothing more. If the former, it's tricky because ideally you'd have more observations per cycle. Can you show me what the Fourier transform looks like? And what the smoothed values from the differentiator look like? Also, about the data: what is "ground truth |
Beta Was this translation helpful? Give feedback.
-
Thanks for your quick answer! Thank you very much for your time and consideration. POINT: all steps below are carried out regarding pre-processed data. stlsq_optimizer = ps.STLSQ(threshold=0.0001, normalize_columns=True)
sr3_optimizer = ps.SR3(threshold=0.1, thresholder="l2")
poly_library = ps.PolynomialLibrary(degree=4, include_interaction=True)
sfd_differentiate = ps.SmoothedFiniteDifference()
model = ps.SINDy(optimizer=stlsq_optimizer,
feature_library=poly_library,
feature_names=['u', 'v'],
differentiation_method=sfd_differentiate)
model.fit(x_train, t=t_train)
model.print(lhs=['du/dt', 'dv/dt'], precision=4)
print('Model score: %f' %
model.score(x_train, t=t_train))
x_sim = model.simulate(x0_train, t_train, integrator='odeint')
x_dot_true = model.differentiate(x_train, t=t_train)
x_dot_predicted = model.predict(x_train)x_dot_true :SmoothedFiniteDifferencenp.array(<AxesArray([[ 1.71223759e+01, -9.27486533e+00],
import scipy as sp
import matplotlib.pyplot as plt
from scipy.fft import fftfreq
from scipy.fft import fft, ifft, fft2, ifft2
# Read your data
df = pd.read_excel(r'C:\Users\istock\Desktop\Projects\ii_resampled.xlsx')
x_train = df.drop(['Time'], axis=1).to_numpy()
t_train = df.Time.values
T = t_train[-1] # seconds
N = x_train.shape[0] # measurements
dt = np.diff(t_train)[0]
f = fftfreq(len(t_train), np.diff(t_train)[0])
u = x_train[:, 0]
v = x_train[:, 1]
u_FFT = fft(u)
v_FFT = fft(v)
# Find the peak frequency
max_idx_u = np.argmax(np.abs(u_FFT))
max_idx_v = np.argmax(np.abs(v_FFT))
peak_freq_u = abs(f[max_idx_u])
peak_freq_v = abs(f[max_idx_v])
print("Peak frequency for u: ", peak_freq)
print("Peak frequency for v: ", peak_freq)
plt.figure(figsize=(6,6), dpi=150)
plt.subplot(2,1,1)
plt.plot(f[:N//2], np.abs(u_FFT[:N//2]), label='U')
plt.ylabel('|$\hat{u}_n$|', fontsize=20)
plt.legend()
plt.subplot(2,1,2)
plt.plot(f[:N//2], np.abs(v_FFT[:N//2]), color='r', label='V')
plt.xlabel('$f_n$ [$s^{-1}$]', fontsize=20)
plt.ylabel('|$\hat{v}_n$|', fontsize=20)
plt.legend()
plt.show()
|
Beta Was this translation helpful? Give feedback.
-
|
I did not get similar plots. When I plot the Fourier transform of the data you've given me, I see a lot more higher frequencies. I re-smoothed it the data to get rid of them. import derivative
import pysindy as ps
import matplotlib.pyplot as plt
import numpy as np
from scipy.signal import blackman
from scipy.fft import fftfreq, fft
data = np.array(<your data>)
tf2 = ps.SmoothedFiniteDifference(smoother_kws = {"window_length":20})
tf2(data[:, 0])
N = len(data)
w = blackman(N)
x_axis = fftfreq(N)[:N//2]
plt.plot(x_axis, np.abs(fft(w * data[:,0])[:N//2]), label="your data")
plt.plot(x_axis, np.abs(fft(w * tf2.smoothed_x_)[:N//2]), label="additional savgol")
plt.legend()
tf1 = derivative.Kalman(alpha=1e2)
tf2 = ps.SmoothedFiniteDifference(smoother_kws = {"window_length":20})
plt.plot(tf1.x(data[:500, 0], t=np.arange(500)), label="kalman")
tf2(data[:500, 0])
plt.plot(tf2.smoothed_x_, label="savgol")
plt.plot(data[:500,0], ".", label="your data")
plt.legend()If the dynamics don't truly include those highest frequencies, try ramping up the power of the SG filter or using something like total variation. That said, if the higher frequencies are real and you're sampling near the nyquist rate of your dynamics, it'll be very hard to recover any true equations. You might also try using the discrete SINDy method in this case (
|
Beta Was this translation helpful? Give feedback.
-
|
Thank you, Jacob! 1- Which data have you used to produce these plots? 2- I ran the codes you sent but : 3- Would you please elaborate on the first plot? Are you comparing the Fourier Transform of the data with the FT of its derivatives? |
Beta Was this translation helpful? Give feedback.
-
|
Unfortunately, Github does not allow me to attach the Original Data and the Pre-processed one inside a spoiler. It gives this error: Original Data
Pre-processed Data
|
Beta Was this translation helpful? Give feedback.
-
Hi there,
I hope this message finds you well. I am writing to inform you that I am working on a thesis related to a turbulent shear layer and have obtained 2D LDA-derived time-series data for this purpose. My research is centered on discovering governing equations using the SINDy algorithm. However, I have encountered a challenge as I have been unable to obtain a satisfactory answer with SINDy despite over a year of attempts.
The original data set comprises 5000 samples with four key features:
1- U_vel [m/s]
2- U_timestamps [ms]
3- V_vel [m/s] and
4- V_timestamps [ms].
Notably, the U & V velocity data is not recorded simultaneously as there is a time gap between them. However, I have successfully addressed this issue using preprocessing techniques such as Resampling and Imputation.
By applying the Pivot concept in Python, I combined the U_vel and V_vel features into a single timestamp column feature. However, SINDy's results were undesirable, with the given equations extremely underfitting. I have tried varying hyperparameters, optimizers, etc to address this issue, but nothing has worked.
I try to include all the necessary information related to my work for your understanding (while attaching the Original Data inside a spoiler, I get the "Unresponsive Error" on Github). I would appreciate your assistance in overcoming the processing challenges.
Best regards,
Reza
Reproducing Code:
Beta Was this translation helpful? Give feedback.
All reactions