Methods

Source data and code

Forecast and observed data were sourced from the European COVID-19 Forecast Hub, available to view at https://covid19forecasthub.eu/ . All Hub data are now archived at:

Data for this work were downloaded on 30th May 2023. These data are available in the Github repository for this paper at: https://github.com/epiforecasts/eval-by-method/tree/main/data

The codebase for this paper is publicly available at:

Comments and code contributions are welcome - please use Github Issues.

Please cite code using:

Study participation

Forecasting teams were recruited to the European Covid-19 Forecast Hub using existing networks and ECDC publicity. Any forecaster was eligible to participate, and there were no selection criteria. All participation was voluntary and unrenumerated. Forecasters contributed a standard set of metadata describing their team and model, and uploaded forecasts weekly. Forecasters were optionally able to express uncertainty by reporting a distribution of up to 23 probabilistic quantiles for each prediction. Forecasts were validated against minimal formatting requirements for quantile intervals and that values were positive integers.

For this study, we collected all forecasts from between 8 March 2021 to 10 March 2023. We excluded forecasts of hospitalisations, which experienced multiple changes in source data during the study period. We excluded forecasts that did not report the full set of 23 quantiles, in order to ensure fair comparison among probabilistic model results. We also excluded baseline and ensemble models created by the Hub team.

Eligibility criteria for models contributing case (left) and death (right) forecasts to the European COVID-19 Forecast Hub, March 2021 - March 2023

Characteristics of contributing models:

Model characteristics contributing to the European COVID-19 Forecast Hub, by method used, number of countries targeted, and number of forecasts contributed.
Model Method Country Targets Case forecasts Death forecasts
AMM-EpiInvert Statistical Multi-country 2,788
CovidMetrics-epiBATS Statistical Single-country 343
DSMPG-bayes Semi-mechanistic Multi-country 760
EuroCOVIDhub-baseline Statistical Multi-country 13,082 13,040
FIAS_FZJ-Epi1Ger Mechanistic Single-country 264 264
GoeWroc-BaseBayes Semi-mechanistic Single-country 12
HZI-AgeExtendedSEIR Mechanistic Single-country 382 382
ICM-agentModel Agent-based Single-country 334 334
IEM_Health-CovidProject Mechanistic Multi-country 7,710 7,708
ILM-EKF Semi-mechanistic Multi-country 11,998 11,961
ITWW-county_repro Semi-mechanistic Multi-country 650 600
Imperial-DeCa Semi-mechanistic Multi-country 571
Imperial-RtI0 Semi-mechanistic Multi-country 571
Imperial-sbkp Semi-mechanistic Multi-country 571
JBUD-HMXK Mechanistic Multi-country 1,324 1,324
KITmetricslab-bivar_branching Statistical Single-country 8
Karlen-pypm Mechanistic Multi-country 3,208 3,186
LANL-GrowthRate Semi-mechanistic Multi-country 3,692 3,696
LeipzigIMISE-SECIR Mechanistic Single-country 16 16
MIMUW-StochSEIR Mechanistic Single-country 76 76
MIT_CovidAnalytics-DELPHI Mechanistic Multi-country 348 500
MOCOS-agent1 Agent-based Single-country 386 386
MUNI-ARIMA Statistical Multi-country 10,979 11,314
MUNI-LaggedRegARIMA Statistical Multi-country 736
MUNI-VAR Statistical Multi-country 976 976
MUNI_DMS-SEIAR Mechanistic Single-country 224 200
PL_GRedlarski-DistrictsSum Mechanistic Single-country 378
RobertWalraven-ESG Statistical Multi-country 9,190 10,465
SDSC_ISG-TrendModel Statistical Multi-country 1,756 1,744
UB-BSLCoV Statistical Single-country 96 96
UC3M-EpiGraph Agent-based Single-country 94
ULZF-SEIRC19SI Mechanistic Single-country 249 249
UMass-MechBayes Mechanistic Multi-country 5,948
UMass-SemiMech Semi-mechanistic Multi-country 1,888 1,904
UNED-PreCoV2 Statistical Single-country 147 147
UNIPV-BayesINGARCHX Statistical Multi-country 426
USC-SIkJalpha Mechanistic Multi-country 12,900 12,688
UpgUmibUsi-MultiBayes Semi-mechanistic Single-country 99 99
bisop-seirfilter Mechanistic Single-country 32 32
bisop-seirfilterlite Mechanistic Multi-country 336 336
epiMOX-SUIHTER Mechanistic Single-country 134 134
epiforecasts-EpiExpert Judgement Multi-country 945 948
epiforecasts-EpiExpert_Rt Judgement Multi-country 404 404
epiforecasts-EpiExpert_direct Judgement Multi-country 394 392
epiforecasts-EpiNow2 Semi-mechanistic Multi-country 8,843 7,721
epiforecasts-weeklygrowth Statistical Multi-country 5,971
itwm-dSEIR Mechanistic Single-country 406 406
prolix-euclidean Semi-mechanistic Multi-country 800 800

Number of models participating in forecasting 1-week ahead case incidence for each country over the study period.

We explored how models selected geographic targets over time. Forecasters both added and removed targets among the set for which they forecast each week. This figure shows number of targets submitted by each model at the 1-week ahead horizon. We noted the same variation for 2-4 week forecasts.

The number of countries targeted by each model forecasting for 1-week ahead incidence of Covid-19 cases or deaths, shown by classification of model structure. Each forecaster selected any number of targets from a set of 32 countries.

Data processing

Model structure classification

Human raters classified models according to methodological structures.

Classification of models by number of raters in total and agreement on model structure
Model Final classification Agreement Total raters Semi-mechanistic Mechanistic Agent-based Statistical Judgement Other Machine Learning
AMM-EpiInvert Statistical FALSE 3 0 0 0 2 0 1 0
DSMPG-bayes Semi-mechanistic FALSE 3 2 0 0 1 0 0 0
ILM-EKF Semi-mechanistic FALSE 3 2 0 0 1 0 0 0
ITWW-county_repro Semi-mechanistic FALSE 3 2 0 0 1 0 0 0
Imperial-DeCa Semi-mechanistic FALSE 3 2 0 0 1 0 0 0
Imperial-RtI0 Semi-mechanistic FALSE 3 2 0 0 1 0 0 0
Imperial-sbkp Semi-mechanistic FALSE 3 2 0 0 1 0 0 0
KITmetricslab-bivar_branching Statistical FALSE 3 1 0 0 2 0 0 0
Karlen-pypm Mechanistic FALSE 3 0 2 0 1 0 0 0
LANL-GrowthRate Semi-mechanistic FALSE 3 2 0 0 1 0 0 0
SDSC_ISG-TrendModel Statistical FALSE 3 0 0 0 2 0 1 0
UMass-SemiMech Semi-mechanistic FALSE 4 3 1 0 0 0 0 0
USC-SIkJalpha Mechanistic FALSE 4 1 3 0 0 0 0 0
UpgUmibUsi-MultiBayes Semi-mechanistic FALSE 3 2 0 0 1 0 0 0
bisop-seirfilter Mechanistic FALSE 4 1 3 0 0 0 0 0
bisop-seirfilterlite Mechanistic FALSE 4 1 3 0 0 0 0 0
prolix-euclidean Semi-mechanistic FALSE 4 3 0 0 0 0 1 0
CovidMetrics-epiBATS Statistical TRUE 3 0 0 0 3 0 0 0
EuroCOVIDhub-baseline Statistical TRUE 4 0 0 0 4 0 0 0
FIAS_FZJ-Epi1Ger Mechanistic TRUE 3 0 3 0 0 0 0 0
GoeWroc-BaseBayes Semi-mechanistic TRUE 3 3 0 0 0 0 0 0
HZI-AgeExtendedSEIR Mechanistic TRUE 3 0 3 0 0 0 0 0
ICM-agentModel Agent-based TRUE 3 0 0 3 0 0 0 0
IEM_Health-CovidProject Mechanistic TRUE 4 0 4 0 0 0 0 0
JBUD-HMXK Mechanistic TRUE 3 0 3 0 0 0 0 0
LeipzigIMISE-SECIR Mechanistic TRUE 3 0 3 0 0 0 0 0
MIMUW-StochSEIR Mechanistic TRUE 3 0 3 0 0 0 0 0
MIT_CovidAnalytics-DELPHI Mechanistic TRUE 3 0 3 0 0 0 0 0
MOCOS-agent1 Agent-based TRUE 3 0 0 3 0 0 0 0
MUNI-ARIMA Statistical TRUE 3 0 0 0 3 0 0 0
MUNI-LaggedRegARIMA Statistical TRUE 3 0 0 0 3 0 0 0
MUNI-VAR Statistical TRUE 3 0 0 0 3 0 0 0
MUNI_DMS-SEIAR Mechanistic TRUE 3 0 3 0 0 0 0 0
PL_GRedlarski-DistrictsSum Mechanistic TRUE 3 0 3 0 0 0 0 0
RobertWalraven-ESG Statistical TRUE 3 0 0 0 3 0 0 0
UB-BSLCoV Statistical TRUE 3 0 0 0 3 0 0 0
UC3M-EpiGraph Agent-based TRUE 3 0 0 3 0 0 0 0
ULZF-SEIRC19SI Mechanistic TRUE 3 0 3 0 0 0 0 0
UMass-MechBayes Mechanistic TRUE 3 0 3 0 0 0 0 0
UNED-PreCoV2 Statistical TRUE 3 0 0 0 3 0 0 0
UNIPV-BayesINGARCHX Statistical TRUE 3 0 0 0 3 0 0 0
epiMOX-SUIHTER Mechanistic TRUE 3 0 3 0 0 0 0 0
epiforecasts-EpiExpert Judgement TRUE 3 0 0 0 0 3 0 0
epiforecasts-EpiExpert_Rt Judgement TRUE 3 0 0 0 0 3 0 0
epiforecasts-EpiExpert_direct Judgement TRUE 3 0 0 0 0 3 0 0
epiforecasts-EpiNow2 Semi-mechanistic TRUE 3 3 0 0 0 0 0 0
epiforecasts-weeklygrowth Statistical TRUE 2 0 0 0 2 0 0 0
itwm-dSEIR Mechanistic TRUE 3 0 3 0 0 0 0 0

Epidemic trend identification

We categorised each week as “Stable”, “Decreasing”, or “Increasing”, based on the difference over a three-week moving average of incidence (with a change of +/-5% as “Stable”).

Trends (cases)

Trends (deaths)

Variant phase identification

Genomic surveillance data were obtained from three sources: ECDC (covering 30 European countries), UKHSA (Great Britain), and the Swiss Federal Office of Public Health (Switzerland). Variant lineages were mapped to six named phases in expected chronological order: Alpha, Delta, Omicron-BA.1, Omicron-BA.2, Omicron-BA.4/5, and Omicron-BQ/XBB. For each country, we identified the first week in which each named variant exceeded 50% of sequenced samples. We enforced chronological ordering by removing any out-of-sequence phases, then expanded phase assignments to all weeks by filling forward and backward from observed transition dates. This per-location approach accounts for the fact that variant dominance dates differed substantially across European countries. Where genomic surveillance data were too sparse to identify a transition (Hungary), we supplemented with epidemiological reports to set the Alpha-to-Delta transition date.

Variant phases identified by dominant variant in each location and week

Model specification

We developed this work in two major versions. The first aimed to evaluate both model structure and geographic target specificity (one or multiple countries) and their effect on model performance. The second substantially revised the aim and scope:

  • Limited aim to developing the use of a modelling approach to evaluate forecast performance
  • Specified a single response (model structure) and treating other covariates as confounding factors
  • Improved range of confounding factors (epidemiological outcomes covering both cases and deaths; adding a covariate for dominant variant phase)
  • More extensive model checking and diagnostics to address highly skewed residuals

See the (Sensitivity analysis){#sensitivity-analysis} below for alternative model specifications.

Covariate selection

We aimed to assess the direct relationship betwen Method (i.e. model structure), and LWIS. We show the assumptions underlying our strategy in a causal diagram. We emphasise that model specification is extremely flexible, and we approach this work is exploratory rather than fully inferential.

Directed acyclic graph of assumed causal relationships among forecast performance (wis), model structure (Method), and covariates.

The joint model therefore adjusts for:

  • Forecast target difficulty (which block backdoor paths from epidemic dynamics into the score), as:
    • EpiTarget, Trend, Horizon, Incidence, Location, VariantPhase
  • Forecaster strategy, as:
    • CountryTargets, Model

Adjustment gives the partial association between model structure and LWIS with these covariates held fixed at their observed values.

Model fitting

We noted our outcome of LWIS had a strongly skewed PDF.

We used a [] distribution [skew, kurtosis].

We used the mgcv package v1.9-4 using R 4.5, with the formula:

wis ~ Epi_target + s(Method, bs = “re”) + s(CountryTargets, bs = “re”) + s(Incidence) + s(Trend, bs = “re”) + s(Location, bs = “re”) + s(VariantPhase, bs = “re”) + s(Horizon, by = Model, k = 3, bs = “sz”) + s(Model, bs = “re”)

We assessed model fit visually plotting observed-versus-fitted values and inspecting residuals in Q-Q and histogram plots.

Results

Effects are deviations from the grand mean for each covariate under a sum-to-zero constraint. Negative values indicate better-than-average performance.

Partial effects across covariates

Partial effect on forecast performance of key covariates included in the fully adjusted model for additional covariates. Partial effects of individual model structure and individual model are presented in the main text.

  1. Spatial covariates

  1. Temporal covariates

Partial effects on the raw scale

The main-text figures and tables report exponentiated effects (multiplicative ratios relative to the grand-mean WIS). The table below gives the underlying raw partial effects from the fully adjusted model.

Raw partial effects on the log-WIS scale from the fully adjusted generalised additive mixed model. Negative values indicate better-than-average performance. 95% CI = 95% confidence interval.
Variable Group Partial effect (95% CI), log scale
CountryTargets Multi-country 0 (-0.002, 0.002)
Single-country 0 (-0.002, 0.002)
Epi_target Deaths -1.84 (-1.867, -1.813)
Location AT -0.041 (-0.107, 0.025)
BE 0.171 (0.109, 0.234)
BG -0.341 (-0.413, -0.269)
CH 0.131 (0.068, 0.193)
CY 0.151 (0.087, 0.215)
CZ 0.003 (-0.062, 0.068)
DE -0.259 (-0.329, -0.189)
DK 0.063 (-0.001, 0.128)
EE -0.112 (-0.178, -0.046)
ES 0.002 (-0.065, 0.068)
FI 0.244 (0.182, 0.305)
FR 0.077 (0.011, 0.144)
GB -0.057 (-0.122, 0.008)
GR 0.091 (0.028, 0.153)
HR -0.056 (-0.123, 0.011)
HU 0.086 (0.023, 0.15)
IE 0.158 (0.096, 0.221)
IS 0.35 (0.289, 0.411)
IT -0.213 (-0.283, -0.143)
LI -0.078 (-0.142, -0.013)
LT -0.163 (-0.231, -0.095)
LU 0.139 (0.077, 0.201)
LV -0.149 (-0.217, -0.081)
MT 0.055 (-0.008, 0.118)
NL 0.191 (0.127, 0.255)
NO -0.23 (-0.3, -0.161)
PL -0.295 (-0.364, -0.227)
PT -0.023 (-0.087, 0.04)
RO -0.016 (-0.081, 0.05)
SE -0.056 (-0.123, 0.011)
SI 0.006 (-0.059, 0.071)
SK 0.172 (0.11, 0.233)
Method Agent-based 0 (-0.001, 0.001)
Judgement 0 (-0.001, 0.001)
Mechanistic 0 (-0.001, 0.001)
Semi-mechanistic 0 (-0.001, 0.001)
Statistical 0 (-0.001, 0.001)
Model AMM-EpiInvert -0.023 (-0.102, 0.056)
CovidMetrics-epiBATS -0.247 (-0.44, -0.055)
DSMPG-bayes -0.192 (-0.35, -0.033)
FIAS_FZJ-Epi1Ger 0.004 (-0.181, 0.189)
GoeWroc-BaseBayes 0.274 (-0.115, 0.663)
HZI-AgeExtendedSEIR -0.148 (-0.325, 0.03)
ICM-agentModel 0.122 (-0.029, 0.274)
IEM_Health-CovidProject -0.238 (-0.317, -0.159)
ILM-EKF 0.048 (-0.027, 0.122)
ITWW-county_repro -0.006 (-0.139, 0.127)
Imperial-DeCa 0 (-0.421, 0.421)
Imperial-RtI0 0 (-0.421, 0.421)
Imperial-sbkp 0 (-0.421, 0.421)
JBUD-HMXK -0.21 (-0.305, -0.115)
KITmetricslab-bivar_branching -0.028 (-0.421, 0.364)
Karlen-pypm -0.187 (-0.274, -0.1)
LANL-GrowthRate -0.033 (-0.116, 0.05)
LeipzigIMISE-SECIR -0.015 (-0.405, 0.374)
MIMUW-StochSEIR 0.385 (0.133, 0.637)
MIT_CovidAnalytics-DELPHI 0.276 (0.119, 0.433)
MOCOS-agent1 -0.286 (-0.475, -0.097)
MUNI-ARIMA -0.189 (-0.264, -0.113)
MUNI-LaggedRegARIMA 0.119 (-0.105, 0.343)
MUNI-VAR 0.214 (0.118, 0.309)
MUNI_DMS-SEIAR -0.237 (-0.42, -0.055)
PL_GRedlarski-DistrictsSum -0.235 (-0.42, -0.049)
RobertWalraven-ESG 0.093 (0.018, 0.168)
SDSC_ISG-TrendModel 0 (-0.421, 0.421)
UB-BSLCoV 0.08 (-0.158, 0.318)
UC3M-EpiGraph 0.002 (-0.304, 0.308)
ULZF-SEIRC19SI -0.318 (-0.49, -0.147)
UMass-MechBayes -0.035 (-0.16, 0.089)
UMass-SemiMech 0.175 (0.087, 0.262)
UNED-PreCoV2 -0.102 (-0.3, 0.097)
UNIPV-BayesINGARCHX 0.326 (0.184, 0.468)
USC-SIkJalpha 0.177 (0.103, 0.251)
UpgUmibUsi-MultiBayes 0 (-0.421, 0.421)
bisop-seirfilter 0.092 (-0.258, 0.441)
bisop-seirfilterlite 0.119 (-0.051, 0.289)
epiMOX-SUIHTER 0.029 (-0.342, 0.401)
epiforecasts-EpiExpert -0.041 (-0.166, 0.084)
epiforecasts-EpiExpert_Rt -0.131 (-0.307, 0.045)
epiforecasts-EpiExpert_direct -0.001 (-0.159, 0.157)
epiforecasts-EpiNow2 0.151 (0.076, 0.227)
epiforecasts-weeklygrowth -0.192 (-0.268, -0.115)
itwm-dSEIR -0.088 (-0.265, 0.088)
prolix-euclidean 0.497 (0.414, 0.579)
Trend Decreasing 0.147 (-0.326, 0.62)
Increasing 0.325 (-0.148, 0.798)
Stable -0.471 (-0.945, 0.002)
VariantPhase Alpha -0.419 (-0.627, -0.211)
Delta -0.15 (-0.357, 0.057)
Omicron-BA.1 0.3 (0.092, 0.507)
Omicron-BA.2 0.203 (-0.004, 0.411)
Omicron-BA.4/5 0.067 (-0.141, 0.274)
Omicron-BQ/XBB 0 (-0.208, 0.207)

Sensitivity analyses

Data processing

Population normalisation

We scored forecasts and observations normalised by country population per 100,000. We also scored against the raw total incident weekly count per target. This had near-zero impact on model fitting or results on the log-transformed scores. Difference in median WIS (across all targets):

Target Scale Raw count Per 100k
Cases log 0.292 0.278
Cases natural 2426 34.4
Deaths log 0.239 0.093
Deaths natural 14.7 0.199

Natural scale WIS

We present results from scoring forecast error on the natural scale (difference between observation and prediction in count of case or death incidence), from which the WIS is then calculated and the analysis repeated.

Characteristics of models and forecasts sampled from the European COVID-19 Forecast Hub, March 2021-2023. Models (%) shows number of models and percentage of all included models. Single-country shows models targeting one country as a fraction and percentage of models in each method group.
Cases
Deaths
Models (%) Models (%) Single-country (%)
Overall 42 (100%) 38 (100%) NA
Method
Semi-mechanistic 9 (21.4%) 10 (26.3%) 2/12 (17%)
Statistical 11 (26.2%) 7 (18.4%) 4/12 (33%)
Mechanistic 16 (38.1%) 16 (42.1%) 10/17 (59%)
Agent-based 3 (7.1%) 2 (5.3%) 3/3 (100%)
Judgement 3 (7.1%) 3 (7.9%) 0/3 (0%)
Geographic scope
Single-country 19 (45.2%) 14 (36.8%) NA
Multi-country 23 (54.8%) 24 (63.2%) NA

Model specification

We tried several alternative model specifications.

Covariate selection

Domain Variable Model version Impact
Forecast target
EpiTarget V2* None
VariantPhase V2* Some
Trend All Untested
Horizon All Untested
Incidence All Untested
Location All Untested
Forecaster strategy
CountryTargets All Untested
Model All Untested
Outcome Method All

* Covariates included in V2 were added based on reviewer feedback.

Fitting

We tried a Gaussian distribution with a log link. Because LWIS is strongly right-skewed, this resulted in high residual skewness (5.5). We avoid using additional data transformations (e.g. an additional log transform) on the outcome (LWIS) as this would violate propriety of the score.

The fitted relationships and partial effects remained stable across parameterisations.