如何解决Python GEKKO-如何在ODE中使用数组中的值
我们有一个项目,我们确实需要帮助。
基本上,我们正在尝试使用GEKKO解决一个多方程组。但是,参数(miu)之一是由神经网络预测的。 但是,当我们尝试将预测的数据和方程式放在一起时,会出现多个错误。
import numpy as np
from gekko import GEKKO,brain
import pandas as pd
import matplotlib.pyplot as plt
from math import e
m = GEKKO(remote=False) # create GEKKO model -- optimization and accesses solvers of constrained,unconstrained,continuous,and discrete problems
KdQ = 0.001 #degree of degradation of glutamine (1/h)
mG = 1.1e-12# 1.1e-10 #glucose maintenance coefficient (mmol/cell/hour)
YAQ = 0.1#0.90 #yield of ammonia from glutamine
YLG = 0.1 #2 #yield of lactate from glucose
YXG = 2.2e8 #yield of cells from glucose (cells/mmol)
YXQ = 0.5e9#1.5e9 #yield of cells from glutamine (cells/mmol)
KL = 150 #lactate saturation constant (mM)
KA = 40 #ammonia saturation constant (mM)
Kdmax = 0.01 #maximum death rate (1/h)
mumax = 0.044 #maximum growth rate (1/h)
KG = 30#1 #glucose saturation constant (mM)
KQ = 0.22 #glutamine saturation constant (mM)
mQ = 0 #glutamine maintenance coefficient (mmol/cell/hour)
kmu = 0.01 #intrinsic death rate (1/h)
Klysis = 2e-2 #rate of cell lysis (1/h)
Ci_star = 100 #inhibitor saturation concentration (mM)
qi = 2.5e-10 #specific inhibitor production rate (1/h)
#Flow,volume and concentration
Fo = 0 #Feed-rate (L/h)
Fi = 0 #Feed-rate (L/h)
V = 3 #volume (L)
SG = 653 #glucose concentration in the Feed (mM)
SQ = 58.8 #glutamine concentration in the feced (mM)
#Load experimental data
from Experimental_Data import tspan,glucose,glutamine,glutamate,lact,ammonia,cell_br1,cell_br2
# create GEKKO parameter
t = np.linspace(0,144,99)
m.time = t
XT= m.Var(value=5e8,name='XT') #total cell density (MMcells/L)
XV = m.Var(value=5e8,lb=0,name='XV') #viable cell density (MMcells/L)
from test_ann import b,x
# mu values are given by neural network
mu2 = b.think(x)
mu1 = np.array(mu2)
#mu = m.abs3(mu2)
mu = m.sos1(mu1)
Kd = m.Intermediate(Kdmax*(kmu/(mu+kmu))) #death rate(1/h)
# create GEEKO equations
m.Equation(XT.dt()== mu*XV )
m.Equation(XV.dt() == ((mu - Kd)*XV ))
# solve ODE
m.options.IMODE = 4 #Simulation #2-Regression mode
m.options.soLVER = 1 #Public software version
m.options.NODES = 3 #Default
m.options.COLDSTART = 2
# objective
m.solve(display=False)
# objective
#m.Obj(sum([ (z[j]-1)**2 + y for j in range(p)]))
#figure,axes = plt.subplots(nrows=5,ncols=1)
plot1 = plt.figure(1)
plt.plot(t,XV.value,label='viable cell')
#axes[0].plot(t,XT.value,label='total cell')
plt.xlabel='Time [hr]'
plt.ylabel='Concentration [cells/ml]'
plt.legend()
plot1 = plt.figure(2)
plt.xlabel='Time [hr]'
plt.ylabel='Concentration [mM]'
plt.legend()
plot1 = plt.figure(3)
plt.plot(tspan,'bx',label = 'Lactate measured')
plt.xlabel='Time [hr]'
plt.ylabel='Concentration [mM]'
plt.legend()
plot1 = plt.figure(4)
plt.plot(tspan,'ro',label = 'Ammonia measured')
plt.plot(tspan,label = 'glutamine measured')
plt.xlabel='Time [hr]'
plt.ylabel='Concentration [mM]'
plt.legend()
plot1 = plt.figure(5)
plt.plot(m.time,mu,label='\u03BC')
plt.plot(m.time,Kd,label='Kd')
plt.xlabel='Time [hr]'
plt.ylabel='Miu[1/h]'
plt.legend()
plt.show()
使用实验数据获取数据
import pandas as pd
#Load experimental data
df = pd.read_excel(r'path')
sheet = df[0:9] #we have to include row 235
tspan = sheet['TIME']
cell_br1= sheet['CELL_BR1']
cell_br2= sheet['CELL_BR2']
由于我无法将excel文件放在此处,因此数据如下:
然后使用此模块(ann_test)预测miu
from gekko import GEKKO
from gekko import brain
import numpy as np
import matplotlib.pyplot as plt
from numpy import diff
from scipy.interpolate import CubicSpline
xm = np.array([ 0.0,23.0,47.0,71.5,95.0,119.0,143.0 ]) # 47.0,deriv1 = 0
from Experimental_Data import cell_br1,cell_br2
def spline(cell):
m = GEKKO()
m.options.IMODE=2
c = [m.FV(value=0) for i in range(4)]
x = m.Param(value=xm)
cell = np.array(cell)
y = m.CV(value=cell)
y.FSTATUS = 1
# polynomial model
m.Equation(y==c[0]+c[1]*x+c[2]*x**2+c[3]*x**3)
c[0].STATUS=1
m.solve(disp=False)
c[1].STATUS=1
m.solve(disp=False)
c[2].STATUS=1
c[3].STATUS=1
m.solve(disp=False)
pbr = [c[3].value[0],c[2].value[0],\
c[1].value[0],c[0].value[0]]
# print(pbr)
xp = np.linspace(0,100)
plot1 = plt.figure(1)
if cell[0] == cell_br2[0]:
plt.plot(xm,cell_br2,'ko',label ='BR2')
plt.plot(xp,np.polyval(pbr,xp),'g:',linewidth=2)
elif cell[0] == cell_br1[0] :
plt.plot(xm,'mo',label ='BR1')
plt.plot(xp,'r:',linewidth=2)
plt.xlabel('time(hr)')
plt.ylabel('cells')
plt.legend()
dx = diff(xp)
dy1 = diff(np.polyval(pbr,xp))
deriv1 = dy1/dx
time =np.linspace(0,99)
plot1 = plt.figure(2)
if cell[0] == cell_br2[0]:
plt.plot(time,deriv1,'b:',linewidth=2,label ='BR2')
elif cell[0] == cell_br1[0]:
plt.plot(time,'m:',label ='BR1')
plt.xlabel('time(hr)')
plt.ylabel('miu(1/h)')
plt.legend()
#plt.show()
return(deriv1)
m = GEKKO()
from Experimental_Data import cell_br1,glucose
b = brain.Brain(remote=True)
b.input_layer(2)
b.layer(linear=5)
b.layer(tanh=3)
b.layer(tanh=5)
b.output_layer(1)
x_s = np.linspace(0,99)
xg = np.array([ 0.0,\
95.0,144.0 ])
cells_spline = CubicSpline(xm,cell_br1)
y_cells = cells_spline(x_s)
miu_1 = spline(cell_br1)
miu_2 = spline(cell_br2)
scale = [1.0e6,1.0e4]
x = (x_s,y_cells) #,y_glucose) #Inputs (3)
y1 = (miu_1) #Output (2)
y2 = (miu_2) #Output (2)
b.learn(x,y1) # train
b.learn(x,y2) # train
yp = b.think(x) # validate
x_1 = np.linspace(0,198)
xp = np.linspace(0,99)
yyp = np.array(yp)
miu = np.reshape(yyp,(99,))
plot1 = plt.figure(3)
plt.plot(x_s,miu,'r-',label = 'Predicted ')
plt.plot(x_s,miu_1,'.',label = 'Experimental points')
plt.xlabel('Time [hr]')
plt.ylabel('miu [1/h]')
plt.legend()
plt.show()
问题是我无法将miu的值(来自an_test)与微分方程合并。
这是我得到的错误:
TypeError:无法根据规则“安全”将数组数据从dtype('O')转换为dtype('float64')
可以请人帮忙吗?
解决方法
问题可能是您正在使用m.sos1()
函数为您的微分方程生成mu
:
mu = m.sos1(mu1)
Kd = m.Intermediate(Kdmax*(kmu/(mu+kmu))) #death rate(1/h)
# create GEEKO equations
m.Equation(XT.dt()== mu*XV )
m.Equation(XV.dt() == ((mu - Kd)*XV ))
要将参数向量(与m.time
相同的长度)转换为微分方程,请使用m.Param()
创建mu
参数。
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。