R: передача параметров в функцию (пакет deSolve)

Question

Недавно я загрузил пакет deSolve для решения ODE-моделей. Я запрограммировал модель в соответствии с некоторым примером кода, который я нашел, но, похоже, проблема связана с передачей параметров, так как я получаю массаж ошибок об неизвестных параметрах.

library(deSolve) 

    model <- function(t, y, parms) { 
    dY1 = -(k1 * y[1]) + (k2 * y[6]) - (k13 * y[1] * 400*sin(((2*pi)/period_ca)*t-phase_ca)+600) + (k14 * y[2]) - (k17 * y[1] * 400*sin(((2*pi)/period_dg)*t-phase_dg)+600) + (k18 *y[11]) - (k3*y[1] * AA) + (k4 * y[7]) 
    dY2 = (k13 * y[1] * 400*sin(((2*pi)/period_ca)*t-phase_ca)+600) - (k14 * y[2]) - (k15 * y[2] *400*sin(((2*pi)/period_dg)*t-phase_dg)+600) + (k16 * y[3]) - (k5 * y[2]) + (k6 * y[8]) - (k7 * y[2] * AA) + (k8 * y[9]) 
    dY3 = (k15 * y[2] * 400*sin(((2*pi)/period_dg)*t-phase_dg)+600) - (k16 * y[3]) - (k9 * y[3]) + (k10 * y[10]) 
    dY4 = -(k11 * y[4]) + (k12 * y[5]) + (k19 * y[11] * AA) - (k20 * y[4]) 
    dY5 = (k11 * y[4]) - (k12 * y[5]) 
    dY6 = (k1 * y[1]) - (k2 * y[6]) 
    dY7 = (k3 + y[1] * AA) - (k4 * y[7]) 
    dY8 = (k5 * y[2]) - (k6 * y[8]) 
    dY9 = (k7 * y[2] * AA) - (k8 * y[9]) 
    dY10 = (k9 * y[3]) - (k10 * y[10]) 
    dY11 = (k17 * y[1] * 400*sin(((2*pi)/period_dg)*t-phase_dg)+600) - (k18 * y[11]) - (k19 * y[11] * AA) + (k20 * y[4]) 
    list(c(dY1, dY2, dY3,dY4, dY5, dY6,dY7, dY8, dY9, dY10, dY11)) 
    } 

    yini <- c(y1 = 1000, y2 = 0, y3 = 0, y4 = 0, y5 = 0, y6 = 20, y7 = 0, y8 = 0, y9 = 0, y10 = 0, y11 = 0) 
    times <- seq(from = 0, to = 5000, by = 0.1) 
    parms <- c(AA=11000, k1=1, k2=50, k3=1.2e-7, k4=0.1, k5=1.2705, k6=3.5026, k7=1.2e-6, k8=0.1, k9=1, k10=0.1, k11=2, k12=0.2, k13=0.0006, k14=0.5, k15=7.998e-6, 
      k16=8.6348, k17=6e-7, k18=0.1, k19=1.8e-5, k20=2, period_ca=100, phase_ca=0, period_dg=100, 
      phase_dg=0) 
    out <- ode (times = times, y = yini, func = model, parms = parms) 

Здесь dY1 - dY11 представляют собой дифференциальные уравнения для некоторых компонентов системы. Parms - вектор, определяющий необходимые значения параметров, yini определяет начальные условия и временную шкалу.

Я получаю следующее сообщение об ошибке:

Ошибка в FUNC (время, состояние, Parms, ...): объект 'k1' не найден

Я совершенно новой для R и не понимают происхождение проблемы (весь код кода, который я нашел, был построен таким же образом).

Answer 1

Вы должны найти несколько примеров из пакета deSolve. Чтобы использовать параметры в функции, вы должны использовать with функции:

model <- function(t, y, parms) { 
    with(as.list(c(y, parms)), { 
    dY1 = ... 
    . 
    . 
    . 
    dY11 = ... 

    list(c(dY1, dY2, dY3,dY4, dY5, dY6,dY7, dY8, dY9, dY10, dY11)) 
}) 
} 

Тем не менее я получил ошибку, потому что интеграция не удается, но это происходит не потому, что код является неправильным. Возможно, это не так на вашем ПК.

Вы можете попробовать этот код для решения ОДУ, что дает некоторые предупреждения, но интеграция успешна (на моем компьютере):

out <- ode(times = times, y = yini, func = model, parms = parms, method = "bdf")