用零替换模型(ODE 系统)中的负值

Question

我目前正在使用 deSolve 求解常微分方程组，并且想知道是否有任何方法可以防止微分变量值低于零。我看过其他一些关于在向量、数据帧等中将负值设置为零的帖子，但由于这是一个生物模型(T 细胞计数变为负值是没有意义的)，我需要从一开始就阻止它发生，这样这些值就不会扭曲结果，而不仅仅是替换最终输出中的负数。

Answer 1

我的标准方法是将状态变量转换为无约束的尺度。对于正变量，最明显/标准的方法是写下 log(x) 的动力学方程。而不是 x .

例如，对于传染病流行的易感-感染-恢复 (SIR) 模型，其中方程为 dS/dt = -beta*S*I; dI/dt = beta*S*I-gamma*I; dR/dt = gamma*I我们会天真地把梯度函数写成

gfun <- function(time, y, params) {
   g <- with(as.list(c(y,params)),
       c(-beta*S*I,
          beta*S*I-gamma*I,
          gamma*I)
       )
   return(list(g))
}

如果我们做log(I)而不是I是状态变量(原则上我们也可以用 S 来做到这一点，但实际上 S 接近边界的可能性要小得多)，那么我们有 d(log(I))/dt = (dI/dt)/I = beta*S-gamma ;其余方程需要使用exp(logI)引用I 。所以:

gfun_log <- function(time, y, params) {
   g <- with(as.list(c(y,params)),
       c(-beta*S*exp(logI),
          beta*S-gamma,
          gamma*exp(logI))
       )
   return(list(g))
}

(计算 exp(logI) 一次并存储/重复使用它比计算两次更有效......)