Python GEKKO - ODE'lerimde dizideki değerleri nasıl kullanabilirim

Nov 13 2020

Bir projemiz var ve gerçekten yardıma ihtiyacımız var.

Temelde yapmaya çalıştığımız şey, GEKKO kullanarak çoklu bir denklem sistemini çözmek. Bununla birlikte, parametrelerden biri (miu) sinir ağları tarafından tahmin edilir. Bununla birlikte, tahmin edilen verileri ve denklemleri bir araya getirmeye çalıştığımızda, birden çok hata elde ederiz.

İki programım var: Bu birincisi, asıl olan:

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,lact,'bx', label = 'Lactate measured')


plt.xlabel='Time [hr]' 
plt.ylabel='Concentration [mM]'
plt.legend()


plot1 = plt.figure(4)

plt.plot(tspan,ammonia,'ro', label = 'Ammonia measured')
plt.plot(tspan,glutamine,'bx', 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()

Veriler, Experimental_Data kullanılarak elde edilir

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 dosyasını buraya koyamadığım için veriler aşağıdaki gibidir:

Ve miu bu modül kullanılarak tahmin edilir (ann_test)

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,144,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,cell_br1, 'mo', label ='BR1')
        plt.plot(xp,np.polyval(pbr,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,144,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,deriv1,'m:',linewidth=2, 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, cell_br2, 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,144,99)
xg = np.array([ 0.0 , 23.0 , 47.0 , 71.5 ,\
                95.0 , 119.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,144,198)
xp = np.linspace(0,144,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()

Sorun şu ki miu değerlerini (ann_test'ten) diferansiyel denklemlerle birleştiremiyorum.

Bu elde ettiğim hata:

TypeError: 'safe' kuralına göre dizi verileri dtype ('O') 'dan dtype (' float64 ')' e dönüştürülemiyor

Lütfen birisi yardım edebilir mi?

Yanıtlar

1 JohnHedengren Nov 26 2020 at 11:00

Sorun, diferansiyel denklemleriniz m.sos1()için oluşturmak üzere işlevi kullanıyor olmanız olabilir 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.timeBir diferansiyel denklem içine bir parametre vektörü (aynı uzunlukta ) m.Param()almak için, muparametreyi oluşturmak için kullanın .