pomp 包 COVID SEIR 模型最小二乘误差回溯

pomp package COVID SEIR model least squares errors traceback

我尝试为英国建立 SEIR 模型以评估实施的遏制措施,并在此处找到了一些带有 pomp 包的代码:https://kingaa.github.io/clim-dis/parest/parest.html 我试图将其转移到我的案例中,该案例增加了一个阶段 (E) 和另外三个变量。最后我想做一个最小二乘估计来找到最佳贝塔。 Data_UK_beta0 由变量日期(从 0 到 165 的整数)和 new_cases(来自约翰霍普金斯大学数据集)组成。

Data_UK_pomp_beta0 <- pomp(
  data= Data_UK_beta0,
  times ="date", t0=0,
  skeleton = vectorfield(
    Csnippet("
      DS=-beta1*S*I/N;
      DE= beta1*S*I/N-delta1*E;
      DI=delta1*E-(gamma1+eta1)*I;
      DR=gamma1*I;")),
  rinit = Csnippet("
      S=S_0;
      E=E_0;
      I=I_0;
      R=N-S_0-E_0-I_0;"),
  statenames = c("S","E","I","R"),
  paramnames = c("beta1","delta1","gamma1","eta1","N","S_0","E_0","I_0")) 


sse_UK_beta0 <- function (params) {
  x <- trajectory(Data_UK_pomp_beta0,params=params)
  discrep <- x["I",,]-obs(Data_UK_pomp_beta0)
  sum(discrep^2)
}


install.packages("apricom")
library(apricom)

beta_reg <- function (beta0) {
  params <- c(beta1=beta0, delta1=1/5.1, gamma1=1, eta1=0.012649, N=67886004, S_0=67886004, E_0=5, I_0=2)
  sse(params)
}


beta0 <- seq(from=1,to=40,by=1)
SSE <- sapply(beta0, beta_reg)

得到以下错误回溯(不幸的是德语,但我想消息应该很清楚:)

  Fehler bei der Auswertung des Argumentes 'x' bei der Methodenauswahl für Funktion 'as.matrix': Argument "dataset" fehlt (ohne Standardwert) 
7.
h(simpleError(msg, call)) 
6.
.handleSimpleError(function (cond) 
.Internal(C_tryCatchHelper(addr, 1L, cond)), "Argument \"dataset\" fehlt (ohne Standardwert)", 
    base::quote(as.matrix(dataset))) 
5.
as.matrix(dataset) 
4.
sse(params) 
3.
FUN(X[[i]], ...) 
2.
lapply(X = X, FUN = FUN, ...) 
1.
sapply(beta0, beta_reg) 

我做错了什么?

您从 apricom 导入的 sse 函数与此问题无关(据我所知)。 (这也与 C(++) 代码编译没有任何关系,因此您问题中的 [compiler-errors] 标记有点误导。)

你没有给我们一个方法来获取你的 Data_UK_beta0 数据集,所以我无法重现这个,但我假设你真的想要这样的东西:

beta_reg <- function (beta0) {
    params <- c(beta1=beta0, delta1=1/5.1, gamma1=1, 
                eta1=0.012649, N=67886004, S_0=67886004, E_0=5, I_0=2)
    sse_UK_beta0(params)
}
beta0 <- seq(from=1,to=40,by=1)
SSE <- sapply(beta0, beta_reg)

您还应该意识到,在您的工作中还有许多其他潜在的棘手问题。一个突然出现的问题是,您可能将模型中的 流行率 感染(I 隔间中的当前数字)与 病例报告 数据,它是 发病率 的(滞后的、有偏差的、不完善的)衡量指标,即每单位时间的新病例数 ...