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

我正在研究使用的常微分方程系统deSolve,并想知道是否有任何方法可以防止差分变量值低于零.我已经看过一些关于在矢量,数据框等中将负值设置为零的其他帖子,但由于这是一个生物模型(并且T细胞计数变为负值没有意义),需要阻止它开始发生,所以这些值不会扭曲结果,而不仅仅是替换最终输出中的底片.

我的标准方法是将状态变量转换为无约束的比例.对正变量执行此操作的最明显/标准方法是写下动力学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)而不是我是状态变量(原则上我们也可以用S做这个,但实际上S不太可能接近边界),那么我们就有了d(log(I))/dt = (dI/dt)/I = beta-gamma; 其余的方程式需要exp(logI)用来指代I.所以:

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

(计算exp(logI)一次并存储/重新使用它而不是计算它两次会稍微高效......)