여러 매개 변수가있는 ODEINT (시간 종속)

Question

ODEINT를 사용하여 단일 일차 주문 ODE를 해결하려고합니다. 다음은 코드입니다. 나는 3 개의 시간 포인트에 대해 y의 3 개의 값을 얻을 것으로 기대한다. 내가 고민하는 문제는 mt 및 nt의 n 번째 값을 dydt으로 계산하는 것입니다. ODEINT는 반복에 따라 mt 및 nt 대신 3 번째 값 (0, 1 또는 2)을 모두 전달한다고 생각합니다. 이 때문에, I이 오류 얻을 : 흥미롭게

RuntimeError: The size of the array returned by func (4) does not match the size of y0 (1).

를 I과 같이 단일 값 인 초기 상태를 대체 (되어야 있으면) : A0 = [2] * 4 코드가 작동하지만 나에게 솔루션으로 4X4 매트릭스를 제공합니다.

mt = np.array([3,7,4,2]) # Array of constants 
nt = np.array([5,1,9,3]) # Array of constants 
c1,c2,c3 = [-0.3,1.4,-0.5] # co-efficients 
para = [mt,nt] # Packing parameters 

#Test ODE function 
def test (y,t,extra): 
    m,n = extra 
    dydt = c1*c2*m - c1*y - c3*n 
    return dydt 

a0= [2] # Initial Condition 
tspan = range(len(mt)) # Define tspan 

#Solving the ODE 
yt= odeint(test, a0,tspan,args=(para,)) 
#Plotting the ODE 
plt.plot(tspan,yt,'g') 
plt.title('Multiple Parameters Test') 
plt.xlabel('Time') 
plt.ylabel('Magnitude') 

첫 번째 미분 방정식은 다음과 같습니다

dy/dt = c1*(c2*mt-y(t)) - c3*nt

이 방정식은 내가 모델링하려고 쥐의 내분비 시스템의 일부를 나타냅니다. 이 시스템은 2 탱크 시스템과 유사합니다. 첫 번째 탱크는 [알 수없는 속도로] 특정 호르몬을 받지만 센서는 특정 시간 간격 (1 초)마다 해당 레벨 (mt)을 감지합니다. 이 탱크는 두 번째 탱크로 공급되며이 호르몬의 레벨은 다른 센서에 의해 감지됩니다. (y). 레벨을 감지하는 센서가 서로 독립적이며 서로 보정하지 않기 때문에 별도의 변수를 사용하여 레벨을 분류했습니다. 'c2'는 두 수준 간의 상관 관계를 나타내는 공효로 간주 될 수 있습니다. 또한,이 호르몬의 탱크 1에서 탱크 2 로의 이동은 확산에 의해 유도됩니다. 이 호르몬은 생화학 적 과정에 의해 더 많이 소모됩니다 (두 번째 탱크의 배수 밸브와 유사 함). 현재 어떤 매개 변수가 소비에 영향을 미치는지는 분명하지 않습니다. 그러나 다른 센서가 특정 시간 간격 (이 경우에도 1 초)에 소비되는 호르몬 양을 감지 할 수 있습니다.

따라서 mt 및 nt은 특정 시점에서 호르몬의 농도/수준입니다. 비록 코드의 길이가 4 요소 일지라도,이 배열은 나의 연구에서 훨씬 길다. 모든 센서는 1 초 간격으로 농도를보고합니다. 따라서 tspan은 1 초 간격으로 시간 점으로 구성됩니다.

목적은 두 번째 탱크 (y)에서이 호르몬의 농도를 수학적으로 결정한 다음 실험 데이터를 기반으로 이러한 계수의 값을 최적화하는 것입니다. 나는이 배열 mt과 nt을 정의 된 ODE에 전달하고 문제없이 MATLAB에서 ODE45를 사용하여 해결할 수있었습니다. 나는이 RunTimeError를 실행하면서 Python으로 코드를 복제하려고 시도했다.

Answer 1

일반적인 미분 방정식을 사용하여이 시스템을 모델링하려면 주석 시간에 m 및 n 사이의 값을 가정해야합니다. 하나의 가능한 모델은 선형 보간법을 사용하는 것입니다. 여기 scipy.interpolate.interp1d을 사용하여 및 nfunc(t)이라는 예제를 작성하고 mt 및 nt 샘플을 기반으로하는 스크립트가 있습니다.여기

import numpy as np 
from scipy.integrate import odeint 
from scipy.interpolate import interp1d 
import matplotlib.pyplot as plt 

mt = np.array([3,7,4,2]) # Array of constants 
nt = np.array([5,1,9,3]) # Array of constants 

c1, c2, c3 = [-0.3, 1.4, -0.5] # co-efficients 

# Create linear interpolators for m(t) and n(t). 
sample_times = np.arange(len(mt)) 
mfunc = interp1d(sample_times, mt, bounds_error=False, fill_value="extrapolate") 
nfunc = interp1d(sample_times, nt, bounds_error=False, fill_value="extrapolate") 

# Test ODE function 
def test (y, t): 
    dydt = c1*c2*mfunc(t) - c1*y - c3*nfunc(t) 
    return dydt 

a0 = [2] # Initial Condition 
tspan = np.linspace(0, sample_times.max(), 8*len(sample_times)+1) 
#tspan = sample_times 

# Solving the ODE 
yt = odeint(test, a0, tspan) 

# Plotting the ODE 
plt.plot(tspan, yt, 'g') 
plt.title('Multiple Parameters Test') 
plt.xlabel('Time') 
plt.ylabel('Magnitude') 
plt.show() 

스크립트에 의해 생성 된 플롯이다

참고 대신 단지 sample_times의 용액을 생성 (즉, 시간 0, 1, 2시, 3), I tspan을 밀도가 높은 점 집합으로 설정하십시오. 이것은 샘플 시간 사이의 모델의 동작을 보여줍니다.