search

cmmid / hhh4-workshop

Workshop on hhh4 and the endemic-epidemic framework
https://cmmid.github.io/hhh4-workshop/
15 stars 4 forks source link

Paste your results for one-week ahead forecasts of rotavirus in Berlin #1

Open jbracher opened 2 years ago

jbracher commented 2 years ago

You can paste the return of colMeans(scores(owa)) here, if you like along with your control argument or a little piece of R code. See my example post below.

jbracher commented 2 years ago

Example:

library("hhh4addon")
data("rotaBE")

ctrl <- list(end = list(f = addSeason2formula(~0 + fe(1, unitSpecific = TRUE),
                                              S = 1)),
             ar = list(f = ~ 0 + fe(1, unitSpecific = TRUE)),
             family = "NegBin1",
             subset = 6:312)

fit <- hhh4(rotaBE, ctrl)
owa <- oneStepAhead(fit, tp = c(312, 363))
colMeans(scores(owa))
     logs       rps       dss       ses 
 1.931189  1.639351  3.139681 22.410975
bastistician commented 2 years ago 
library("hhh4addon")
data("rotaBE")

## let the power law act on (o+1) => AR is integrated in NE
rotaBE1 <- rotaBE
neighbourhood(rotaBE1) <- neighbourhood(rotaBE1) + 1

form1 <- ~ 0 + fe(1, unitSpecific = TRUE)
ctrl <- list(end = list(f = addSeason2formula(form1, S = 2)),
             ne = list(f = addSeason2formula(form1, S = 2),
                       weights = W_powerlaw(maxlag = max(neighbourhood(rotaBE1)))),
             family = "NegBinM", subset = 6:312)
fit <- hhh4(rotaBE1, ctrl)
plot(fit, total = TRUE)
owa <- oneStepAhead(fit, tp = c(312, 363))
colMeans(scores(owa))
     logs       rps       dss       ses 
 1.879683  1.587602  2.932538 21.733318

Did anyone achieve something smaller?

I observed that some more complicated models do not converge at all time points during the rolling forecasts. One could use oneStepAhead(..., type="first") in that case, i.e. produce all forecasts from the initial fit at the beginning of the test period without re-fitting.

jbracher commented 2 years ago 
library("hhh4addon")
data("rotaBE")

## let the power law act on (o+1) => AR is integrated in NE
rotaBE1 <- rotaBE
neighbourhood(rotaBE1) <- neighbourhood(rotaBE1) + 1

form1 <- ~ 0 + fe(1, unitSpecific = TRUE)
ctrl <- list(end = list(f = addSeason2formula(form1, S = 2)),
             ne = list(f = addSeason2formula(form1, S = 2),
                       weights = W_powerlaw(maxlag = max(neighbourhood(rotaBE1)))),
             family = "NegBinM", subset = 6:312)
fit <- profile_par_lag(rotaBE1, ctrl)
plot(fit, total = TRUE)
owa <- oneStepAhead_hhh4lag(fit, tp = c(312, 363))
colMeans(scores(owa))
#        logs       rps       dss       ses 
# 1.864334  1.550379  2.702720 20.628763

Just a little tweak ;)