Covid-19 增长率（Bootstrapping/Time 系列）
Covid-19 Growth rate (Bootstrapping/Time Series)
我正在尝试编写 R 代码以获得 COVID-19 的增长率。
可以在插入的图像上找到方程式，其中 i(t) 是时间 t 时的感染人数。
我想我可以编写代码，例如，如果等式是简单的增长率 i(t) = i0 * exp(r*t)；然而，我所指的论文实现了以下方法：“每个时间点的时间序列数据，使用 bootstrap 提取 1000 次重复观察的方法，假设它遵循泊松分布，其中 ^Data 为平均值."
我以前从未编写过这样的方程式:(我希望这是成为更好的编码员和统计学家的绝佳机会和启蒙。
有人可以帮我解决这个问题吗？我真的很感激任何类型的提示或建议。我从不期待完整的代码或任何类似的东西；我很想找到一个可以很好地开始这件事的资源。
这是一个使用一些艾滋病病例数据的示例（摘自 Dobson 和 Barnett 的 GLM 介绍）。
aids_dat <- read.csv("https://raw.githubusercontent.com/bbolker/goettingen_2019/master/data/aids.csv")
## generate more useful time variable
aids_dat <- transform(aids_dat, time = year - min(year) + (quarter-1)/4,
cumcases = cumsum(cases))
nobs <- nrow(aids_dat)
log-linear适合数据
（你可以在这里使用任何你想要的模型）
f1 <- lm(log(cases) ~ time, data = aids_dat)
生成模拟数据集
set.seed(101)
fitted <- predict(f1)
s <- replicate(1000,
rpois(nobs, lambda = exp(fitted)),
simplify = FALSE)
这从 log-linear 拟合中获取预测值并添加泊松变化。乔威尔等人。用我不明白的累积病例数做一堆事情（而且通常认为这不是一个好主意）；看起来他们生成预测值，创建累积案例计数，并立即将其差值以返回案例？值得小心这些东西（记住流行和发病率之间的区别），但我不太确定 Chowell 等人是什么。正在做。
改装模型
refit <- lapply(s,
function(x) lm(log(x+0.1) ~ time, data = aids_dat))
我使用 log(x+0.1) 以防泊松样本生成零值。 （最好用Poisson GLM之类的。）
获取每个改装模型的预测值
boot_ensemble <- sapply(refit, predict)
您还可以使用 *apply 或 for 循环从每个改装模型中提取系数，计算流行病指标（增长率或 R0）等——无论您想要什么值 bootstrap 分布为.
情节
png("boot.png")
plot(cases ~ time, data = aids_dat)
matlines(aids_dat$time, exp(boot_ensemble),
type = "l", lty = 1,
col = adjustcolor("black", alpha = 0.01))
matlines(aids_dat$time,
col = "red", lty = 1,
t(apply(exp(boot_ensemble), 1, range)))
dev.off()
您可以通过将 quantile 应用于预测集合的每一行来获得逐点置信区间，但请参阅 Juul 等人 2021 年的警告说明。 （您可以使用 fda::fbplot() 函数绘制与他们相似的曲线。）
最后的警告，我怀疑使用泊松（而不是例如负二项式）参数 bootstrapping 通常会低估不确定性。
Juul、Jonas L.、Kaare Græsbøll、Lasse Engbo Christiansen 和 Sune Lehmann。 “Fixed-Time 描述性统计低估了流行曲线整体的极端情况。”自然物理 17，没有。 1（2021 年 1 月）：5-8。 https://doi.org/10.1038/s41567-020-01121-y.
我正在尝试编写 R 代码以获得 COVID-19 的增长率。
可以在插入的图像上找到方程式，其中 i(t) 是时间 t 时的感染人数。
我想我可以编写代码，例如，如果等式是简单的增长率 i(t) = i0 * exp(r*t)；然而，我所指的论文实现了以下方法：“每个时间点的时间序列数据，使用 bootstrap 提取 1000 次重复观察的方法，假设它遵循泊松分布，其中 ^Data 为平均值."
我以前从未编写过这样的方程式:(我希望这是成为更好的编码员和统计学家的绝佳机会和启蒙。
有人可以帮我解决这个问题吗？我真的很感激任何类型的提示或建议。我从不期待完整的代码或任何类似的东西；我很想找到一个可以很好地开始这件事的资源。
这是一个使用一些艾滋病病例数据的示例（摘自 Dobson 和 Barnett 的 GLM 介绍）。
aids_dat <- read.csv("https://raw.githubusercontent.com/bbolker/goettingen_2019/master/data/aids.csv")
## generate more useful time variable
aids_dat <- transform(aids_dat, time = year - min(year) + (quarter-1)/4,
cumcases = cumsum(cases))
nobs <- nrow(aids_dat)
log-linear适合数据
（你可以在这里使用任何你想要的模型）
f1 <- lm(log(cases) ~ time, data = aids_dat)
生成模拟数据集
set.seed(101)
fitted <- predict(f1)
s <- replicate(1000,
rpois(nobs, lambda = exp(fitted)),
simplify = FALSE)
这从 log-linear 拟合中获取预测值并添加泊松变化。乔威尔等人。用我不明白的累积病例数做一堆事情（而且通常认为这不是一个好主意）；看起来他们生成预测值，创建累积案例计数，并立即将其差值以返回案例？值得小心这些东西（记住流行和发病率之间的区别），但我不太确定 Chowell 等人是什么。正在做。
改装模型
refit <- lapply(s,
function(x) lm(log(x+0.1) ~ time, data = aids_dat))
我使用 log(x+0.1) 以防泊松样本生成零值。 （最好用Poisson GLM之类的。）
获取每个改装模型的预测值
boot_ensemble <- sapply(refit, predict)
您还可以使用 *apply 或 for 循环从每个改装模型中提取系数，计算流行病指标（增长率或 R0）等——无论您想要什么值 bootstrap 分布为.
情节
png("boot.png")
plot(cases ~ time, data = aids_dat)
matlines(aids_dat$time, exp(boot_ensemble),
type = "l", lty = 1,
col = adjustcolor("black", alpha = 0.01))
matlines(aids_dat$time,
col = "red", lty = 1,
t(apply(exp(boot_ensemble), 1, range)))
dev.off()
您可以通过将 quantile 应用于预测集合的每一行来获得逐点置信区间，但请参阅 Juul 等人 2021 年的警告说明。 （您可以使用 fda::fbplot() 函数绘制与他们相似的曲线。）
最后的警告，我怀疑使用泊松（而不是例如负二项式）参数 bootstrapping 通常会低估不确定性。
Juul、Jonas L.、Kaare Græsbøll、Lasse Engbo Christiansen 和 Sune Lehmann。 “Fixed-Time 描述性统计低估了流行曲线整体的极端情况。”自然物理 17，没有。 1（2021 年 1 月）：5-8。 https://doi.org/10.1038/s41567-020-01121-y.