使用 python GEKKO 进行预测

Question

我已经使用 sysid 函数创建了一个 ARX 模型，现在我想将该模型推断为其他输入。

我首先创建了这个模型并确定了其训练数据子集的参数。

m = GEKKO(remote=False)
na = n_a
nb = n_b
yp_train, p, K = m.sysid(t_train, u_train, y_train, na, nb)

让这个模型适合我有什么方法可以使用相同的参数并利用不同的输入来使用它们吗？

这就是我想要的

def arxouriço2(df, inputs, outputs, n_a, n_b, train_slice=None, test_slice=None):
    datafreime_arx = df.copy()
    if isinstance(inputs, pd.Series):
        inputs = pd.DataFrame(inputs)
    if isinstance(outputs, pd.Series):
        outputs = pd.DataFrame(outputs)
    colunas_arx = list(inputs.columns) + list(outputs.columns)
    mean_col = []
    std_dev_col = []

    for i in colunas_arx:
        mean = np.nanmean(datafreime_arx[i])
        mean_col.append(mean)
        datafreime_arx[i] -= mean

        std_dev = np.nanstd(datafreime_arx[i])
        std_dev_col.append(std_dev)
        datafreime_arx[i] /= std_dev

    if train_slice is not None:
        train_data = datafreime_arx.iloc[train_slice]
        t_train = train_data['Tempo']
        u_train = inputs.iloc[train_slice]
        y_train = outputs.iloc[train_slice]
    else:
        t_train = datafreime_arx['Tempo']
        u_train = inputs
        y_train = outputs

    if test_slice is not None:
        test_data = datafreime_arx.iloc[test_slice]
        t_test = test_data['Tempo']
        u_test = inputs.iloc[test_slice]
        y_test = outputs.iloc[test_slice]
    else:
        t_test = datafreime_arx['Tempo']
        u_test = inputs
        y_test = outputs

    m = GEKKO(remote=False)
    na = n_a
    nb = n_b
    yp_train, p, K = m.sysid(t_train, u_train, y_train, na, nb)

    for i, col in enumerate(colunas_arx):
        datafreime_arx[col] /= std_dev_col[i]
        datafreime_arx[col] -= mean_col[i]

    resis = y_train - yp_train
    resis2 = resis ** 2
    print('Sum of squared residuals: ' + str(np.sum(resis2)))

    if test_slice is not None:
        yp_test = m.sysid(t_test, u_test, y_test, na, nb)[0]
        plt.figure(figsize=(8, 5))
        plt.subplot(2, 1, 1)
        plt.plot(t_train, u_train)
        plt.legend(inputs.columns)
        plt.ylabel('Temperature (ºC)')
        plt.subplot(2, 1, 2)
        plt.plot(t_test, y_test, label=outputs.columns)
        plt.plot(t_test, yp_test, label=outputs.columns + ' predicted')
        plt.legend()
        plt.ylabel('Rotation (º)')
        plt.xlabel('Time')
        plt.tight_layout()
        plt.show()

    return p, K, m

arxouriço2(df_c3_avg5_c2,df_c3_avg5_c2[['STH3_avg','STH-04']],df_c3_avg5_c2['C14-Y'],10,10,range(465, 829),range(829, 903))

Answer 1

使用

m.arx()

函数创建预测模型。

m.sysid()

函数用于确定ARX模型的参数

p

。

m.arx()

函数使用

p

参数来创建预测。这是一个带有数据的示例：

将 ARX 与

sysid

配合

from gekko import GEKKO
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# load data and parse into columns
url = 'http://apmonitor.com/pdc/uploads/Main/'
data = pd.read_csv(url+'tclab_data3.txt')
t = data['Time']
u = data['Q1']
y = data['T1']

# generate time-series model
m = GEKKO(remote=False)

# system identification
na = 2 # output coefficients
nb = 2 # input coefficients
yp,p,K = m.sysid(t,u,y,na,nb,pred='meas')

plt.figure(figsize=(10,6))
plt.subplot(2,1,1)
plt.plot(t,u)
plt.legend([r'$Q_1$ (%)'])
plt.ylabel('MV Heater (%)')
plt.subplot(2,1,2)
plt.plot(t,y,'b-',label=r'$T_{1,meas}$')
plt.plot(t,yp,'r--',label=r'$T_{1,pred}$')
plt.legend(); plt.ylabel('CV Temp (°C)')
plt.xlabel('Time (sec)')

使用

arx

函数对 ARX 模型进行预测

# load test data
data = pd.read_csv(url+'tclab_data4.txt')
t = data['Time']
um = data['Q1']
ym = data['T1']

# generate time-series model
m = GEKKO(remote=False)

# Build GEKKO ARX model
y,u = m.arx(p)
y[0].value = ym[0]

# load inputs
m.time = t
u[0].value = um

# options
m.options.imode = 4
m.options.nodes = 2

# simulate
m.solve(disp=False)

plt.figure(figsize=(10,6))
plt.subplot(2,1,1)
plt.plot(t,u[0])
plt.legend([r'$Q_1$ (%)'])
plt.ylabel('MV Heater (%)')
plt.subplot(2,1,2)
plt.plot(t,ym,'b-',label=r'$T_{1,meas}$')
plt.plot(t,y[0],'r--',label=r'$T_{1,pred}$')
plt.legend(); plt.ylabel('CV Temp (°C)')
plt.xlabel('Time (sec)')

根据您的功能，修改这部分代码：

if test_slice is not None:
    yp_test = m.sysid(t_test, u_test, y_test, na, nb)[0]

对于这样的事情：

if test_slice is not None:
    mt = GEKKO(remote=False)
    y,u = mt.arx(p)
    mt.time = t_test
    y[0].value = y_test[0] # initial condition
    u[0].value = u_test
    mt.options.IMODE = 4
    mt.options.NODES = 2
    mt.solve(disp=False)
    yp_test = y[0].value