From c8e35ad8bef6686f505ce6675f82d4751dc5e34b Mon Sep 17 00:00:00 2001 From: "blaise.calmel" <blaise.calmel@gmail.com> Date: Fri, 24 Nov 2023 19:57:32 +0100 Subject: [PATCH] [TEST] Test other functions for temporal interpolation --- UI/main.py | 136 ++++++++++++++++++++++++++++++++++++++++------------- 1 file changed, 104 insertions(+), 32 deletions(-) diff --git a/UI/main.py b/UI/main.py index 0583a6e..1ea1900 100644 --- a/UI/main.py +++ b/UI/main.py @@ -13,7 +13,7 @@ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Affero General Public License for more details. You should have received a copy of the GNU Affero General Public License -along with this program. If not, see <https://www.gnu.org/licenses/>. +along with this program. If not, see <https://www.gnu.org/licenses/>. """ import os @@ -26,7 +26,7 @@ import pickle import pandas as pd from matplotlib.backends.backend_qt5agg import NavigationToolbar2QT as NavigationToolbar -from matplotlib.backends import backend_pdf, backend_svg, backend_ps # avoid bug after compilation +# from matplotlib.backends import backend_pdf, backend_svg, backend_ps # avoid bug after compilation import numpy as np from datetime import datetime @@ -1222,7 +1222,7 @@ class ApplicationWindow(QMainWindow, Ui_MainWindow): self.sorted_name = list(self.mean_meas.index) self.columns_checked = columns_checked_new self.sticky_settings.set('columns_checked', self.columns_checked) - self.run_interface() + self.table_meas() def run_interface(self): """ Plot graphs and tables @@ -1240,43 +1240,42 @@ class ApplicationWindow(QMainWindow, Ui_MainWindow): self.freeze = False print('complete') + """essaie avec une fonction souple en DQ(t) (genre un polynôme ou une spline ou un fit non + paramétrique genre LOWESS?), sans x au début, puis on verra pour ajouter une propagation selon x""" + + # Define values + import numpy as np df = self.transects_df.loc[:, ['tr_q_total', 'meas_label_user', 'start_time', 'end_time']] df['mid_time'] = (df['start_time'] + df['end_time']) / 2 df['mid_time'] = df['mid_time'] - np.nanmin(df['mid_time']) df['meas_label_user'] = pd.to_numeric(df['meas_label_user']) mean_q = np.nanmean(df['tr_q_total']) df['delta_q'] = mean_q - df['tr_q_total'] - from sklearn import linear_model - x = df[['meas_label_user', 'mid_time']] - y = df['delta_q'] - regr = linear_model.LinearRegression() - regr.fit(x, y) - print('Intercept: \n', regr.intercept_) - print('Coefficients: \n', regr.coef_) + X = df[['meas_label_user', 'mid_time']].values.reshape(-1, 2) + DQ = df['delta_q'] + x = X[:, 0] + t = X[:, 1] + dq = DQ + + # x = X[:, 0] + # y = X[:, 1] + # z = Y - import statsmodels.api as sm - x = sm.add_constant(x) # adding a constant - model = sm.OLS(y, x).fit() - predictions = model.predict(x) - print_model = model.summary() - print(print_model) + ############################## 3D LINEAR REGRESSION ############################## import matplotlib.pyplot as plt from sklearn import linear_model from mpl_toolkits.mplot3d import Axes3D - X = df[['meas_label_user', 'mid_time']].values.reshape(-1, 2) - Y = df['delta_q'] - x = X[:, 0] - y = X[:, 1] - z = Y + + x_pred = np.linspace(min(x), max(x), 30) # range of porosity values - y_pred = np.linspace(min(y), max(y), 30) # range of brittleness values - xx_pred, yy_pred = np.meshgrid(x_pred, y_pred) - model_viz = np.array([xx_pred.flatten(), yy_pred.flatten()]).T + t_pred = np.linspace(min(t), max(t), 30) # range of brittleness values + xx_pred, tt_pred = np.meshgrid(x_pred, t_pred) + model_viz = np.array([xx_pred.flatten(), tt_pred.flatten()]).T ols = linear_model.LinearRegression() - model = ols.fit(X, Y) + model = ols.fit(X, DQ) predicted = model.predict(model_viz) - r2 = model.score(X, Y) + r2 = model.score(X, DQ) plt.style.use('default') fig = plt.figure(figsize=(12, 4)) @@ -1285,8 +1284,8 @@ class ApplicationWindow(QMainWindow, Ui_MainWindow): ax3 = fig.add_subplot(133, projection='3d') axes = [ax1, ax2, ax3] for ax in axes: - ax.plot(x, y, z, color='k', zorder=15, linestyle='none', marker='o', alpha=0.1) - ax.scatter(xx_pred.flatten(), yy_pred.flatten(), predicted, facecolor=(0, 0, 0, 0), s=20, + ax.plot(x, t, dq, color='k', zorder=15, linestyle='none', marker='o', alpha=0.1) + ax.scatter(xx_pred.flatten(), tt_pred.flatten(), predicted, facecolor=(0, 0, 0, 0), s=20, edgecolor='#70b3f0') ax.set_xlabel('Distance (m)', fontsize=12) ax.set_ylabel('Time (s)', fontsize=12) @@ -1298,8 +1297,81 @@ class ApplicationWindow(QMainWindow, Ui_MainWindow): ax3.view_init(elev=60, azim=165) fig.suptitle('$R^2 = %.2f$' % r2, fontsize=20) fig.tight_layout() - fig.savefig(r"D:\0_INRAE\regression.png", dpi=300, bbox_inches='tight') + fig.savefig(r"D:\0_INRAE\linear_regression.png", dpi=300, bbox_inches='tight') + + ############################## CURVE FIT INTERPOLATION ############################## + from scipy.optimize.minpack import curve_fit + fit_funcs = {"linear": lambda x, a: a * x ** self.exponent, "power": lambda x, a, b: a * x ** b, } + fit_func = 'power' + + popt, pcov = curve_fit(fit_funcs[fit_func], t, dq) + ############################## POLY FIT INTERPOLATION ############################## + eq = np.polyfit(t, dq, 3) + f = np.poly1d(eq) + t_new = np.linspace(min(t), max(t), 50) + q_new = f(t_new) + + fig = plt.figure(figsize=(12, 4)) + ax = fig.add_subplot() + ax.scatter(t, dq, facecolor=(0, 0, 0, 0), s=20, edgecolor='#70b3f0') + ax.plot(t_new, q_new, color='k', zorder=15, alpha=0.9) + ax.set_xlabel('Time (s)', fontsize=12) + ax.set_ylabel('Delta Q (m3/s)', fontsize=12) + fig.savefig(r"D:\0_INRAE\polyfit_regression.png", dpi=300, bbox_inches='tight') + + ############################## SLINEAR / QUADRATIC / CUBIC FUNCTIONS ############################## + + import numpy as np + from scipy.interpolate import interp1d + import matplotlib.pyplot as plt + # Interpolation for different methods: + interpolations_methods = ['slinear', 'quadratic', 'cubic'] + alpha = np.linspace(min(t), max(t), 1000) + + interpolated_points = {} + for method in interpolations_methods: + interpolator = interp1d(t, dq, kind=method, axis=0) + interpolated_points[method] = interpolator(alpha) + + fig = plt.figure(figsize=(12, 4)) + ax = fig.add_subplot() + ax.scatter(t, dq, facecolor=(0, 0, 0, 0), s=20, edgecolor='#70b3f0', label='original points') + for method_name, curve in interpolated_points.items(): + ax.plot(alpha, curve, '-', label=method_name) + # ax.plot(t_new, q_new, color='k', zorder=15, alpha=0.9) + ax.set_xlabel('Time (s)', fontsize=12) + ax.set_ylabel('Delta Q (m3/s)', fontsize=12) + ax.legend() + fig.savefig(r"D:\0_INRAE\slinear_quadratic_cubic__regression.png", dpi=300, bbox_inches='tight') + + ############################## LOWESS FUNCTION ############################## + + import matplotlib.pyplot as plt + from scipy.interpolate import interp1d + import statsmodels.api as sm + # lowess will return our "smoothed" data with a y value for at every x-value + lowess = sm.nonparametric.lowess(dq, t, frac=.3) + # unpack the lowess smoothed points to their values + lowess_t = list(zip(*lowess))[0] + lowess_dq = list(zip(*lowess))[1] + # run scipy's interpolation + f = interp1d(lowess_t, lowess_dq, bounds_error=False) + t_new = np.linspace(min(t), max(t), 1000) + # this this generate y values for our xvalues by our interpolator + # it will MISS values outsite of the x window + # There might be a better approach, but you can run a for loop + # and if the value is out of the range, use f(min(lowess_x)) or f(max(lowess_x)) + dq_new = f(t_new) + fig = plt.figure(figsize=(12, 4)) + ax = fig.add_subplot() + ax.scatter(t, dq, facecolor=(0, 0, 0, 0), s=20, edgecolor='#70b3f0', label='original points') + ax.plot(lowess_t, lowess_dq, '*', color='r', zorder=10, alpha=0.9, label='lowess points') + ax.plot(t_new, dq_new, color='k', zorder=15, alpha=0.9, label='lowess interpolation') + ax.set_xlabel('Time (s)', fontsize=12) + ax.set_ylabel('Delta Q (m3/s)', fontsize=12) + ax.legend() + fig.savefig(r"D:\0_INRAE\lowess_regression.png", dpi=300, bbox_inches='tight') def load_tables(self): """ Draw tables. @@ -1882,9 +1954,9 @@ class ApplicationWindow(QMainWindow, Ui_MainWindow): tbl.setHorizontalHeaderLabels(header_title) header = tbl.horizontalHeader() - for i in range(ncols-1): - header.setSectionResizeMode(i, QtWidgets.QHeaderView.ResizeToContents) - header.setSectionResizeMode(ncols-1, QtWidgets.QHeaderView.Stretch) + header.setSectionResizeMode(0, QtWidgets.QHeaderView.Stretch) + for i in range(1, ncols): + header.setSectionResizeMode(i, QtWidgets.QHeaderView.Interactive) tbl.setVerticalHeaderLabels(self.mean_meas.index) tbl.setSelectionBehavior(QtWidgets.QTableView.SelectRows) -- GitLab