|
13 | 13 | log.basicConfig(level=log.INFO) |
14 | 14 | logging = log.getLogger(__name__) |
15 | 15 |
|
| 16 | +try: |
| 17 | + from statsmodels.tsa.ar_model import AR |
| 18 | + from sklearn.metrics import mean_squared_error |
| 19 | +except (ImportError, SystemError): |
| 20 | + AR = None |
| 21 | + logging.info('statsmodels and sklearn is not installed. Cannot use "rtoc.estimate".') |
16 | 22 |
|
17 | 23 | def lsfit(self, x, y, func="linear", x0=None, n=None): |
18 | 24 | """ |
@@ -75,7 +81,7 @@ def resample(self, x, y, n=None): |
75 | 81 | x (list): x-values (a linspace) |
76 | 82 | y (list): Interpolated y-values |
77 | 83 | """ |
78 | | - if n is None: |
| 84 | + if n is None and self is not None: |
79 | 85 | n = self.config['global']['recordLength'] |
80 | 86 | if len(x) == len(y): |
81 | 87 | x2 = np.linspace(x[0], x[-1], n) |
@@ -290,3 +296,128 @@ def norm(x, y, max=1, min=0, oldMin=None, oldMax=None): |
290 | 296 | y = (y-oldMin)/oldMax |
291 | 297 | y = y*max+min |
292 | 298 | return x, y |
| 299 | + |
| 300 | +def predict(coef, history): |
| 301 | + """ |
| 302 | + Make a prediction give regression coefficients and lag obs |
| 303 | + """ |
| 304 | + yhat = coef[0] |
| 305 | + for i in range(1, len(coef)): |
| 306 | + yhat += coef[i] * history[-i] |
| 307 | + return yhat |
| 308 | + |
| 309 | +def estimate(x,y, n): |
| 310 | + """ |
| 311 | + Estimate n values in future for a signal using ARMA. You need to install the following python packages with pip3: ``pip3 install statsmodels scikit-learn scikit-metrics patsy`` |
| 312 | +
|
| 313 | + Args: |
| 314 | + x (list): List of x-values |
| 315 | + y (list): List of y-values |
| 316 | + n (int): Number of values to be estimated |
| 317 | +
|
| 318 | + Returns: |
| 319 | + x (list): List of estimated x-values |
| 320 | + y (list): List of estimated y-values |
| 321 | + """ |
| 322 | + if AR is None: |
| 323 | + return |
| 324 | + |
| 325 | + x, y = resample(None, x, y, n=len(x)) |
| 326 | + samplerate = 1/(x[1]-x[0]) |
| 327 | + |
| 328 | + #train = np.diff(y) |
| 329 | + train = np.array(y) |
| 330 | + test = np.linspace(x[-1]+1/samplerate, x[-1]+(1/samplerate)*n, n) |
| 331 | + # size = int(len(data) * 0.66) |
| 332 | + # train, test = data[0:size], data[size:] |
| 333 | + # train autoregression |
| 334 | + model = AR(train) |
| 335 | + model_fit = model.fit(maxlag=6, disp=False) |
| 336 | + window = model_fit.k_ar |
| 337 | + coef = model_fit.params |
| 338 | + # walk forward over time steps in test |
| 339 | + history = [train[i] for i in range(len(train))] |
| 340 | + predictions = list() |
| 341 | + for t in range(len(test)): |
| 342 | + yhat = predict(coef, history) |
| 343 | + predictions.append(yhat) |
| 344 | + history.append(yhat) |
| 345 | + error = mean_squared_error(test, predictions) |
| 346 | + print('Test MSE: %.3f' % error) |
| 347 | + |
| 348 | + return test, predictions |
| 349 | + |
| 350 | +def correlate(x1,y1,x2,y2, mode='valid'): |
| 351 | + """ |
| 352 | + Crosscorrelate two signals using numpy.correlate |
| 353 | +
|
| 354 | + Args: |
| 355 | + x1 (list): List of x-values of signal 1 |
| 356 | + y1 (list): List of y-values of signal 1 |
| 357 | + x1 (list): List of x-values of signal 2 |
| 358 | + y1 (list): List of y-values of signal 2 |
| 359 | + mode: {'valid', 'same', 'full') Refer to the convolve docstring. Note that the default is ‘valid’, unlike convolve, which uses ‘full’. |
| 360 | +
|
| 361 | + Returns: |
| 362 | + x (list): List of x-values |
| 363 | + y (list): List of correlated y-values |
| 364 | + """ |
| 365 | + maxL = max([len(y1), len(y2)]) |
| 366 | + x,ys = combine(None, [[x1,y1],[x2,y2]], n=maxL) |
| 367 | + y = np.correlate(ys[0], ys[1], mode) |
| 368 | + |
| 369 | + return x, y |
| 370 | + |
| 371 | +def forwardEuler(x,y, f, U_0, samplerate, T): |
| 372 | + dt = 1/samplerate |
| 373 | + N_t = int(round(float(T)/dt)) |
| 374 | + u = np.zeros(N_t+1) |
| 375 | + t = np.linspace(0, N_t*dt, len(u)) |
| 376 | + u[0] = U_0 |
| 377 | + for n in range(N_t): |
| 378 | + u[n+1] = u[n] + dt*f(u[n], t[n]) |
| 379 | + return t,u |
| 380 | + |
| 381 | +def CCN(): |
| 382 | + """ |
| 383 | + https://machinelearningmastery.com/how-to-get-started-with-deep-learning-for-time-series-forecasting-7-day-mini-course/ |
| 384 | +
|
| 385 | + https://machinelearningmastery.com/time-series-prediction-lstm-recurrent-neural-networks-python-keras/ |
| 386 | +
|
| 387 | + |
| 388 | + Estimate n values in future for a signal using Convolutional Neural Network model. You need to install the following python packages with pip3: ``pip3 install tensorflow keras`` |
| 389 | +
|
| 390 | + Args: |
| 391 | + x (list): List of x-values |
| 392 | + y (list): List of y-values |
| 393 | + n (int): Number of values to be estimated |
| 394 | +
|
| 395 | + Returns: |
| 396 | + x (list): List of estimated x-values |
| 397 | + y (list): List of estimated y-values |
| 398 | + """ |
| 399 | + from keras.models import Sequential |
| 400 | + from keras.layers import Dense |
| 401 | + from keras.layers import Flatten |
| 402 | + from keras.layers.convolutional import Conv1D |
| 403 | + from keras.layers.convolutional import MaxPooling1D |
| 404 | + # define dataset |
| 405 | + X = np.array([[10, 20, 30], [20, 30, 40], [30, 40, 50], [40, 50, 60]]) |
| 406 | + y = np.array([40, 50, 60, 70]) |
| 407 | + # reshape from [samples, timesteps] into [samples, timesteps, features] |
| 408 | + X = X.reshape((X.shape[0], X.shape[1], 1)) |
| 409 | + # define model |
| 410 | + model = Sequential() |
| 411 | + model.add(Conv1D(filters=64, kernel_size=2, activation='relu', input_shape=(3, 1))) |
| 412 | + model.add(MaxPooling1D(pool_size=2)) |
| 413 | + model.add(Flatten()) |
| 414 | + model.add(Dense(50, activation='relu')) |
| 415 | + model.add(Dense(1)) |
| 416 | + model.compile(optimizer='adam', loss='mse') |
| 417 | + # fit model |
| 418 | + model.fit(X, y, epochs=1000, verbose=0) |
| 419 | + # demonstrate prediction |
| 420 | + x_input = np.array([50, 60, 70]) |
| 421 | + x_input = x_input.reshape((1, 3, 1)) |
| 422 | + yhat = model.predict(x_input, verbose=0) |
| 423 | + print(yhat) |
0 commit comments