在 R 中运行 ODE 模型的困难,包括随时间变化的参数(强制函数) [英] Difficulty running an ODE model in R including a parameter that varies by time (forcing function)
问题描述
我正在尝试使用 deSolve 拟合一个相当基本的 ODE 模型,在模型中包括一个随时间变化的参数(感染力;FOI).虽然在没有此参数的情况下运行模型工作正常,但包含时间相关参数会产生错误(见下文).
我对 R 和数学建模比较陌生,并且已经尝试解决这个问题一段时间了.
我已将 FOI 参数创建为值矩阵,然后使用 approxfun 函数进行插值(正如我所见,这适用于强制函数,例如 (v0.2.1) 于 2019 年 3 月 27 日创建上>
I am trying to fit a reasonably basic ODE model using deSolve, including within the model a parameter that varies by time (force of infection; FOI). While running the model without this parameter works fine, inclusion of the time dependent parameter gives an error (see below).
I am relatively new to R and mathematical modelling, and have been trying to solve this problem for some time now.
I have created the FOI parameter as a matrix of values and then used the approxfun function for interpolation (as I have seen this works with forcing functions, e.g. https://rdrr.io/rforge/deSolve/man/forcings.html).
The model without this time-dependent parameter runs without any errors, but attempting to include it gives the error:
Error in checkFunc(Func2, times, y, rho) :
The number of derivatives returned by func() (200) must equal the
length of the initial conditions vector (2)
I cannot figure out how to fix this error, since I only have 2 initial conditions and it seems that including this time-dependent FOI parameter generates many more derivatives.
I am aware that others have asked a similar question, but I have not found this question posed with respect to forcing functions.
Many thanks in advance for any advice.
# Forcing function data
foi <- matrix(ncol=2,byrow=TRUE,data=c(
0, 0.003, 2, 0.03, 3, 0.08, 4,0.1, 5, 0.12, 6, 0.15,
8, 0.16, 10, 0.14,12, 0.12,14,0.08,15, 0.06,16, 0.03,
17, 0.01,18,0.003,19,0.003,20,0.003,30,0.003,40,0.003,
50,0.003,60,0.003,65,0.01, 70,0.08, 72,0.095,74,0.10,
76,0.1, 78,0.08, 80,0.06))
age <- seq(0, 80, by = 1)
input <- approxfun(x = foi[,1], y = foi[,2], method = "constant", rule = 2)
# Function
ab <- function(time, state, pars) {
with(as.list(c(state, pars)), {
import<-c(input(t))
diggP<- (import *iggN) - iggR*iggP
diggN<- (-import*iggN) + iggR*iggP
return(list(c(diggP, diggN)))
})
}
# Initial values
yini <- c(iggP=0, iggN=1)
# Parameters
pars <- c(iggR = 0, import)
# ODE solver
results<- ode(y=yini, times=age, func=foi_model, pars)
I am hoping to make a model in which at each point in time (or in this case age), FOI varies according to the values I have inputted in the FOI matrix. Thus I would like to see how changing FOI over age influences the output of the differential equations.
Your main problem was that you were passing an argument t
to input
, but that variable doesn't exist in your code. Time is passed to your model as an argument called time
. (Also, your model is called ab
not foi_model
, as is stated in the call to ode
, plus pars
doesn't need import
and should be passed to ode
.)
# Load library
library(deSolve)
# Create FOI matrix
foi <- matrix(ncol=2,byrow=TRUE,data=c(
0, 0.003, 2, 0.03, 3, 0.08, 4,0.1, 5, 0.12, 6, 0.15,
8, 0.16, 10, 0.14,12, 0.12,14,0.08,15, 0.06,16, 0.03,
17, 0.01,18,0.003,19,0.003,20,0.003,30,0.003,40,0.003,
50,0.003,60,0.003,65,0.01, 70,0.08, 72,0.095,74,0.10,
76,0.1, 78,0.08, 80,0.06))
# Times for model solution
age <- seq(0, 80, by = 1)
# Linear interpolation function from FOI data
input <- approxfun(x = foi[,1], y = foi[,2], method = "constant", rule = 2)
# Model to be integrated
ab <- function(time, state, parms) {
with(as.list(c(state, parms)), {
##### IMPORTANT #####
import<-input(time) #<- 'time' was previously 't'
#####################
# Derivatives
diggP<- (import *iggN) - iggR*iggP
diggN<- (-import*iggN) + iggR*iggP
# Return results
return(list(c(diggP, diggN)))
})
}
# Initial values
yini <- c(iggP=0, iggN=1)
# Parameters
pars <- c(iggR = 0)
# Solve model
results<- ode(y=yini, times=age, func=ab, parms = pars)
# Plot results
plot(results)
Created on 2019-03-27 by the reprex package (v0.2.1)
这篇关于在 R 中运行 ODE 模型的困难,包括随时间变化的参数(强制函数)的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!