Подбираем нелинейную модель с помощью python и gekko

У меня есть набор деревьев. В этом наборе данных у меня есть уникальный номер графика, последовательность в том порядке, в котором были взяты данные «Измерение», а также среднее значение высоты в метрах и среднего возраста в годах для деревьев. Примерно так: заголовок данных

Затем я определяю модель для прогнозирования роста, используя возраст следующим образом:

Рост = B0 * ((1 - exp (-B1 * Возраст)) ** B2)

Моя цель - определить значения B0, B1 и B2 соответственно. Для этого я использую пакет gekko, чтобы найти параметры моделей со следующим кодом:

num_p = data_gek.Plot.unique()
nmp = 5
number_p = (data_gek.Plot == num_p[nmp])

m = GEKKO()

xm = np.array(data_gek[number_p]['Age'])
x = m.Param(value=xm)

B0 = m.FV(value=38.2) #value=38.2
B0.STATUS = 1

B1 = m.FV(value=0.1) #value=0.1
B1.STATUS = 1

B2 = m.FV(value=2.08) #value=2.08
B2.STATUS = 1

ym =  np.array(data_gek[number_p]['Height'])
z = m.CV(value=ym)

y = m.Var()
m.Equation(y==B0 * ((1 - m.exp(-B1 *x))**B2))
m.Obj(((y-z)/z)**2)

m.options.IMODE = 2
m.options.SOLVER = 1
m.solve(disp=False)

print(B0.value[0],B1.value[0],B2.value[0])
#output 
27.787958561 0.0052435491089 0.21178326158

Однако я не уверен, что поступаю правильно. Можно ли это сделать без начальных значений в параметрах? Потому что я использовал предыдущие значения для B0, B1 и B2 из литературы.

Если вы хотите увидеть мой набор данных и мой процесс, вы можете получить доступ к этой записной книжке в Google Colab < / а>.


python gekko
person Juan Pablo Cuevas    schedule 17.11.2019    source источник

Ответы (1)


arrow_upward
2
arrow_downward

У вашего сценария только одна проблема. Определение z должно быть типа Param или MV вместо CV как z = m.Param(value=ym), потому что это вход для вашей целевой функции.

Вы также можете использовать встроенную цель, если вы определяете y как CV вместо Var. Вам просто нужно включить статус обратной связи FSTATUS=1, чтобы использовать целевую функцию, которая сводит к минимуму разницу между измерениями и прогнозами модели. Вот модифицированная версия вашего скрипта.

регрессионное соответствие

from gekko import GEKKO
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
url = 'http://apmonitor.com/pdc/uploads/Main/data_2nd_order.txt'
data = pd.read_csv(url)
m  = GEKKO()
xm = np.array(data['time'])
x  = m.Param(value=xm)
B0 = m.FV(1); B1 = m.FV(1); B2 = m.FV(1)
B0.STATUS = 1; B1.STATUS = 1; B2.STATUS = 1
ym =  np.array(data['output (y)'])
y  = m.CV(value=ym)
y.FSTATUS = 1
yi = m.Intermediate(B0 * ((1 - m.exp(-B1 *x))**B2))
m.Equation(y==yi)
m.options.IMODE = 2
m.options.SOLVER = 1
m.solve(disp=True)
print(B0.value[0],B1.value[0],B2.value[0])
plt.plot(xm,ym,'ro')
plt.plot(xm,y.value,'b--')
plt.show()
person Eric Hedengren    schedule 18.11.2019
comment
Я пытаюсь использовать ваш код, но он все равно не работает. Я получаю Exception: @error: Solution Not Found. Даже я очищаю свои данные, ища слои и отсутствующие данные. Не могли бы вы взглянуть на [мои данные] (colab.research.google.com / drive /)? Думаю, это сложнее вашего примера. - person Juan Pablo Cuevas; 18.11.2019
comment
Попробуйте изменить начальные значения для B0, B1 и B2 с 1 на начальные значения, которые вы нашли в литературе. B0 = m.FV(38.2) и так далее ... - person reyPanda; 18.11.2019
comment
Еще можно попробовать установить границы для ваших параметров, особенно B1. Большое отрицательное значение приведет к тому, что член m.exp() будет иметь очень большое значение и затруднит схождение решателя. Вы можете добавить границы с помощью B0 = m.FV(lb=0,ub=100). @reyPanda также предлагает хорошие начальные значения. Вы можете комбинировать эти предложения с B0 = m.FV(value=38.2,lb=0,ub=100) - person John Hedengren; 19.11.2019