library(tidyverse); library(knitr); library(kableExtra); library(broom); library(cowplot); 
library(rstan); library(tidybayes); library(scales); library(DT); library(plotly)



theme_set(theme_classic(base_size = 12) + 
            background_grid(color.major = "grey90", 
                            color.minor = "grey95", 
                            minor = "xy", major = "xy") +
            theme(legend.position = "none"))
m <- read_rds("Hierarchical_Model_NegBin.rds")
d <- read_csv("Interactive Model Checking/NegBin/stan_data.csv")
countries <- d %>% distinct(country, country_id)
which_iceland <- d %>% filter(country == "Iceland") %>% .$country_id %>% unique
n_countries <- max(d$country_id)
results <- tidyMCMC(m, conf.int = T, rhat = T, ess = T, 
                    estimate.method = "median", conf.method = "quantile") %>% 
  mutate(par = str_match(term, "[a-zA-Z_2]+")) %>% 
  group_by(par) %>% 
  mutate(num = row_number() %>% as.numeric)

Updates

2020-05-06

  • Conveniently we can define a parameter \(\nabla t\) as the time it takes a country to go from 10% of its maximum diagnosed cases to 90% of its diagnosed cases, and calculate it as \(\nabla t = \frac{\ln(81)}{\beta}\). The country specific time interval and mean time interval posterior distributions are now shown along with the original parameters.

2020-04-14

2020-04-02

  • Model the \(\beta_i\) as lognormal since they are strictly positive.
  • There is enough data to loosen the informative priors, and we have done so.
  • Added posterior predictions for effective reproduction number, \(R_e\), along with some text about it.
  • Reparametrised to use a noncentered parametrisation to aid Stan in performing NUTS.

2020-03-31

  • Reparametrised the prior for S in terms of \(\mu\) the mean and \(\kappa\) the sample size.

2020-03-28

  • Weak priors on \(a_S\) and \(b_S\) to help the NUTS sampler.
  • Rewrote priors for scale parameters.

2020-03-27

  • We now place informative priors on \(\mu_\alpha\), \(\mu_\beta\) and \(\lambda_\phi\) to ensure that the model can be identified.

2020-03-26

  • Implemented a negative binomial likelihood with country-specific overdispersion.
  • Changed filtering criteria for inclusion in model data in order to include more countries.
  • Did some visual inspection of data in order to make sure we don’t lose countries to inner_join strategies. Saw that we were missing USA among others due to different names in different datasets, but is now fixed.

Methods

Parametrisation

Let \(E_i\) and \(I_{i, t}\) be the population and number of infected in country \(i\) at time \(t\). Then the percent of infected can be calculated as

\[ P_{i, t} = \frac{I_{i, t}}{E_i}. \]

In the usual Logistic Regression GLM we could model the percent of infected, as a function of time (in days), with

\[ \log\left(\frac{P_{i, t}}{1 - P_{i, t}}\right) = \alpha_i + \beta_i \cdot t, \]

where \(\alpha_i\) is a measure of how many have been infected in country \(i\) at time \(t = 0\) and \(\beta_i\) is a measure of growth. In the case of COVID-19 infections we don’t know the maximum percent of populations that will be infected, so we have another unknown parameter, the saturation percent at which a country will reach its maximum number of infected, \(S_i\). Thus our model looks like

\[ \log\left(\frac{P_{i, t}}{S_i - P_{i, t}}\right) = \alpha_i + \beta_i \cdot t. \]

These parameters are hard to estimate when data from only one country are used. However, if we were to pool information about them between countries, as in a hierarchical Bayesian model, estimation might be possible. Let

\[ z_{i, t} = \alpha_i + \beta_i \cdot t, \]

where \(\alpha_i\) is a measure of how many have been infected in country \(i\) at time \(t = 0\) and \(\beta_i\) is a measure of growth. Then

\[ P_{i, t} = \frac{S_i}{1 + \exp(-z_{i, t})}, \]

and conditional on some sampling distribution, \(f\), we could write

\[ I_{i, t} \sim \mathrm{f_\theta}(P_{i, t}, E_i), \]

where \(\theta\) contains all relevant parameters to be estimated.

Bayesian Inference

Bayesian inference is a great tool when small amounts of data are to be shared from various sources. In this case the sources are different countries, and the data are cumulative numbers of cases. If we utilize a Negative Binomial likelihood for the observed cumulative cases, then

\[ I_{i, t} \sim \mathrm{NegBin}(P_{i, t} \cdot E_i, \phi_i), \]

where \(\phi_i\) is a country-level effect specifying the amount of overdispersion in that country.

Modeling daily counts

However, there is a lot of implicit correlation in the values of \(I_{i, t}\) for different values of \(t\). Thus, a better parametrisation would be to model the daily number of cases.

Let

\[ z_{i, t} = \alpha_i + \beta_i \cdot t, \]

so that the percent of infected, \(P_i\), is

\[ P_{i, t} = \frac{S_i}{1 + \exp(-z_{i,t})}. \]

If we furthermore write

\[ z^*_{i, t - 1} = \alpha_i + \beta_i \cdot (t - 1), \]

and

\[ P^*_{i, t - 1} = \frac{S_i}{1 + \exp(-z^*_{i, t-1})}, \]

the change in rates between days is

\[ C_{i, t} = P_{i, t} - P^*_{i, t - 1}. \]

Since \(C_{i, t}\) is simply the first derivative of \(P_{i, t}\) with respect to \(t\), we can skip the differencing step and directly model the derivative

\[ C_{i, t} = \frac{d}{dt}P_{i, t} = \beta_i S_i \frac{\exp{(-z_{i, t})}}{(\exp(-z_{i, t}) + 1))^2} \]

Then, conditional on a Negative Binomial likelihood and population size \(E_i\), the daily number of observed cases, \(D_{i, t}\), can be written as

\[ D_{i, t} \sim \mathrm{NegBin}(C_{i, t} \cdot E_i, \phi_i) \]

Parameters

The parameters, \(\alpha\) and \(\beta\), are treated as in a generalized linear model, except that we model \(\beta_i\) on the \(\log\) scale to impose the constraints \(\beta_i > 0\). We put hierarchical priors on \(\alpha_i\) and \(\beta_i\) so that for each country, \(i\),

\[ \begin{aligned} \beta_i &\sim \mathrm{LogNormal}(\mu_\beta, \sigma^2_\beta) \\ \alpha_i &\sim \mathrm{Normal}(\mu_\alpha, \sigma^2_\alpha) \end{aligned} \]

The \(\mu\) parameters are given vaguely informative prior distribution based on our previous Poisson model

\[ \begin{aligned} \mu_\alpha &\sim \mathrm{Normal}(-2.5 , 3^2) \\ \mu_\beta &\sim \mathrm{Normal}(-3, 1^2) \end{aligned} \]

We chose to put priors on the standard deviations, so that \(\sigma_\beta\) and \(\sigma_\alpha\) are given \(\mathrm{Exponential}\) prior distributions

\[ \begin{aligned} \sigma_\alpha &\sim \mathrm{Exponential}(1) \\ \sigma_\beta &\sim \mathrm{Exponential}(2). \end{aligned} \]

The \(S_i\) parameters take on values in \((0, 1)\), so we thought it proper to model them with Beta distributions. By putting hierarchical priors on them we could also share information between countries on the estimated saturation points. Thus the saturation parameter for country \(i\) is sampled as

\[ S_i \sim \mathrm{Beta}(a_S, b_S), \]

where we parametrise in terms of the mean, \(\mu_S\), and prior sample size, \(\kappa_S\)

\[ \begin{aligned} \mu_S &\sim \mathrm{Beta}(1, 99) \\ \kappa_S &\sim \mathrm{Exponential}(0.001). \end{aligned} \]

We can then back-transform into \(a_S\) and \(b_S\) as

\[ \begin{aligned} a_S &= \mu_S \cdot \kappa_S \\ b_S &= (1 - \mu_S) \cdot \kappa_S \end{aligned} \]

We parametrise the negative binomial likelihood in the form of mean and overdispersion. If we write

\[ D_{i, t} \sim \mathrm{NegBin}(\mu_{i, t}, \phi_i), \]

where

\[ \mu_{i, t} = C_{i, t} \cdot E_i, \]

we can write

\[ E[D_{i, t}] = \mu_{i, t} \qquad Var[D_{i, t}] = \mu_{i, t} + \phi_i \cdot \mu_{i, t}^2. \]

Based on Dan Simpson’s excellent post, which is also linked in Stan-Dev’s post on prior choices we put a hierarchical exponential prior on the \(\phi_i\) parameters as follows:

\[ \begin{aligned} \sqrt{\phi_i} &\sim \mathrm{Exponential}(\sigma_\phi) \\ \sigma_\phi &\sim \mathrm{Exponential}(1) \end{aligned} \]

Putting it all together we get

\[ \begin{aligned} D_{i, t} &\sim \mathrm{NegBin}(C_{i, t} \cdot E_i, \phi_i) \\ C_{i, t} = \frac{d}{dt}P_{i, t} &= \beta_i S_i \frac{\exp{(-z_{i, t})}}{(\exp(-z_{i, t}) + 1))^2} \\ z_{i, t} &= \alpha_i + \beta_i \cdot t \\ \beta_i &\sim \mathrm{LogNormal}(\mu_\beta, \sigma^2_\beta) \\ \alpha_i &\sim \mathrm{Normal}(\mu_\alpha, \sigma^2_\alpha) \\ \mu_\beta &\sim \mathrm{Normal}(-3, 1^2) \\ \mu_\alpha &\sim \mathrm{Normal}(-2.5, 3^2) \\ p(\sigma_\alpha) &\sim \mathrm{Exponential}(1) \\ p(\sigma_\beta) &\sim \mathrm{Exponential}(1) \\ \sqrt\phi_i &\sim \mathrm{Exponential}(\sigma_\phi^2)\\ \sigma_\phi &\sim \mathrm{Exponential}(1) \\ S_i &\sim \mathrm{Beta}(a_S, b_S) \\ a_S &= \mu_S \cdot \kappa_S \\ b_S &= (1 - \mu_S) \cdot \kappa_S \\ \mu_S &\sim \mathrm{Beta}(1, 99) \\ \kappa_S &\sim \mathrm{Exponential}(0.001) \end{aligned} \]

Data

We use Icelandic data from the Directorate of Health and supplement with data on other countries from the European CDC, obtained daily.

Data used for modeling were filtered according to

  • Only use data for countries with a population above 200.000.
  • Only keep data for which cases per 1000 inhabitants \(\geq\) 0.1.
  • Only keep data for countries for which the first observation after this filtering is \(\leq\) 0.2
  • Only keep data for countries with more than 6 days of follow-up after the above steps.
  • Countries that are not included for other reasons
    • China
    • South Korea
d %>% 
  group_by(country) %>% 
  summarise(First = min(date),
            Days_In_Data = n(),
            Start_Rate = min(case_rate),
            End_Rate = max(case_rate)) %>% 
  set_names(c("Country", "Entry", "Days in data", "Start", "End")) %>% 
  kable(caption = "Table 1. Summary information about countries used in modeling",
        align = c("l", rep("c", ncol(.) - 1)),
        digits = 3) %>% 
  kable_styling(bootstrap_options = c("striped", "hover")) %>% 
  row_spec(which_iceland, bold = T) %>% 
  add_header_above(c("", "", "", "Rate per 1000" = 2)) %>%
  scroll_box(height = "500px")
Table 1. Summary information about countries used in modeling
Rate per 1000
Country Entry Days in data Start End
Albania 2020-04-04 33 0.106 0.285
Armenia 2020-03-27 41 0.111 0.885
Australia 2020-03-26 42 0.111 0.272
Austria 2020-03-17 51 0.113 1.740
Azerbaijan 2020-04-12 25 0.105 0.205
Bahamas 2020-04-09 28 0.103 0.229
Bahrain 2020-03-14 54 0.128 2.267
Barbados 2020-03-30 38 0.115 0.286
Belarus 2020-04-09 28 0.113 1.941
Belgium 2020-03-18 50 0.108 4.377
Bosnia And Herzegovina 2020-03-31 37 0.107 0.590
Brazil 2020-04-13 24 0.105 0.544
Bulgaria 2020-04-15 22 0.102 0.241
Canada 2020-03-27 41 0.107 1.658
Cape Verde 2020-04-16 21 0.100 0.338
Chile 2020-03-29 39 0.101 1.162
Colombia 2020-04-26 11 0.102 0.171
Costa Rica 2020-04-10 27 0.107 0.150
Croatia 2020-03-26 42 0.101 0.511
Cuba 2020-04-22 15 0.100 0.149
Cyprus 2020-03-25 43 0.105 0.744
Denmark 2020-03-13 55 0.117 1.702
Djibouti 2020-04-08 29 0.124 1.150
Dominican Republic 2020-04-01 36 0.103 0.790
Ecuador 2020-03-29 39 0.106 1.835
Equatorial Guinea 2020-04-25 12 0.156 0.324
Estonia 2020-03-16 52 0.129 1.291
Finland 2020-03-23 45 0.113 0.978
France 2020-03-17 51 0.102 2.042
French Polynesia 2020-03-27 41 0.107 0.215
Gabon 2020-04-29 8 0.110 0.183
Georgia 2020-04-21 16 0.101 0.151
Germany 2020-03-20 48 0.169 1.974
Greece 2020-03-29 39 0.101 0.252
Guinea 2020-04-30 7 0.106 0.142
Hungary 2020-04-09 28 0.101 0.321
Iceland 2020-03-05 62 0.109 5.306
Iran 2020-03-12 56 0.109 1.206
Ireland 2020-03-20 48 0.114 4.502
Israel 2020-03-22 46 0.104 1.912
Italy 2020-03-09 59 0.122 3.518
Jamaica 2020-04-26 11 0.103 0.160
Japan 2020-04-25 12 0.102 0.121
Kazakhstan 2020-04-21 16 0.105 0.231
Kuwait 2020-04-05 32 0.114 1.380
Kyrgyzstan 2020-04-24 13 0.102 0.136
Latvia 2020-03-25 43 0.103 0.470
Lebanon 2020-04-24 13 0.100 0.108
Lithuania 2020-03-27 40 0.108 0.525
Luxembourg 2020-03-16 52 0.125 6.237
Malaysia 2020-04-04 33 0.104 0.200
Maldives 2020-04-21 16 0.130 1.079
Malta 2020-03-19 49 0.109 1.095
Mauritius 2020-03-31 37 0.107 0.277
Mexico 2020-04-25 12 0.101 0.204
Moldova 2020-04-02 35 0.105 1.079
Montenegro 2020-03-27 41 0.107 0.516
Morocco 2020-04-25 12 0.103 0.143
Netherlands 2020-03-19 49 0.120 2.403
New Zealand 2020-03-30 38 0.115 0.238
North Macedonia 2020-03-28 40 0.105 0.732
Norway 2020-03-13 55 0.115 1.469
Oman 2020-04-12 25 0.110 0.550
Palestine 2020-04-30 7 0.102 0.108
Panama 2020-03-25 43 0.104 1.772
Peru 2020-04-09 28 0.134 1.575
Poland 2020-04-06 31 0.108 0.381
Portugal 2020-03-22 45 0.125 2.513
Puerto Rico 2020-04-03 34 0.108 0.656
Qatar 2020-03-14 54 0.113 6.053
Romania 2020-03-31 37 0.101 0.715
Russia 2020-04-13 24 0.108 1.065
Saudi Arabia 2020-04-11 26 0.107 0.883
Serbia 2020-04-01 36 0.103 1.103
Singapore 2020-03-27 41 0.102 3.344
Slovakia 2020-04-08 29 0.106 0.260
Slovenia 2020-03-16 52 0.105 0.695
Spain 2020-03-15 51 0.123 4.693
Sweden 2020-03-16 52 0.103 2.313
Switzerland 2020-03-14 54 0.130 3.483
Turkey 2020-03-30 38 0.110 1.552
Ukraine 2020-04-17 20 0.101 0.289
United Arab Emirates 2020-04-03 34 0.105 1.555
United Kingdom 2020-03-25 43 0.120 2.887
United States 2020-03-23 45 0.107 3.660
Uruguay 2020-04-03 34 0.104 0.194 
p <- d %>% 
  ggplot(aes(days, case_rate, group = country, col = country == "Iceland")) +
  geom_line() +
  scale_y_log10() +
  scale_colour_manual(values = c("grey", "blue")) +
  labs(x = "Days since rate reached 0.02 per 1000",
       y = "Cases per 1000",
       title = "Observed trends for countries used in data",
       subtitle = "Shown as days since a country entered the modeling data")
ggplotly(p)

Software

The model is fit using Stan’s R interface, and the model code can be found in the appendix.

All code is available at https://github.com/bgautijonsson/covid19. We’re sorry for the mess and we’re working on it. We are also going to translate the README files into English.

Results

Convergence

Sampler configuration

  • Four chains were run for 2000 iterations after a warm-up period of 2000 iterations, a total of 4000 iterations. The slowest chain converged in 388 seconds.
  • Adapt_delta was set to 0.99

Warnings

  • There were 7 divergent transitions after warmup.
  • There were 2 transitions after warmup that exceeded maximum treedepth.

We have examined the model carefully to find what causes these divergent transitions, and our current belief is country-specific correlations between \(\beta_i\) and \(S_i\) when the model predicts a high probability of the country being close to its maximum. We have now removed South Korea from the data as the model has a hard time fitting the daily number of new cases in a country after it has reached its asymptote. We are working on steps to overcome this flaw.

Country Level Effects

Conveniently we can define a parameter \(\nabla t\) as the time it takes a country to go from 10% of its maximum diagnosed cases to 90% of its diagnosed cases, and calculate it as \(\nabla t = \frac{\ln(81)}{\beta}\). The country specific posterior distributions are shown below along with the original parameters.

results %>% 
  ungroup %>% 
  filter(par %in% c("beta", "alpha", "S", "phi_inv", "nu", "delta_t")) %>%
  mutate(par = str_replace(par, "phi_inv", "phi")) %>% 
  inner_join(countries, by = c("num" = "country_id")) %>% 
  mutate(par = str_to_title(par)) %>% 
  select(par, country, num, everything(), -num, -term, -std.error) %>% 
  set_names(c("Parameter", "Country", "Median", "Lower", "Upper", "Rhat", "ESS")) %>% 
  kable(digits = 4, align = c("l", "l", rep("c", ncol(.) - 2)),
        caption = "Table 2. Summary of posterior samples of country level effects") %>% 
  kable_styling(bootstrap_options = c("striped", "hover")) %>% 
  add_header_above(c("", "", "", "95% PI" = 2, "Convergence" = 2)) %>% 
  column_spec(1, bold = T) %>% 
  row_spec(which_iceland + c(0, 1, 2) * n_countries, bold = T) %>% 
  collapse_rows(1, valign = "top") %>% 
  scroll_box(height = "600px")
Table 2. Summary of posterior samples of country level effects
95% PI
Convergence
Parameter Country Median Lower Upper Rhat ESS
S Albania 0.0005 0.0003 0.0025 1.0005 2258
Armenia 0.0028 0.0015 0.0073 1.0012 4613
Australia 0.0019 0.0006 0.0058 0.9999 6047
Austria 0.0024 0.0020 0.0035 1.0005 4384
Azerbaijan 0.0006 0.0003 0.0026 1.0018 2168
Bahamas 0.0006 0.0002 0.0036 1.0013 3361
Bahrain 0.0043 0.0028 0.0086 1.0010 5892
Barbados 0.0011 0.0003 0.0048 1.0001 4972
Belarus 0.0042 0.0024 0.0091 0.9999 4997
Belgium 0.0048 0.0044 0.0053 1.0000 15326
Bosnia And Herzegovina 0.0014 0.0008 0.0036 1.0000 2809
Brazil 0.0021 0.0010 0.0064 1.0009 4508
Bulgaria 0.0005 0.0002 0.0033 1.0002 3062
Canada 0.0032 0.0026 0.0049 1.0005 3961
Cape Verde 0.0016 0.0006 0.0055 1.0003 5351
Chile 0.0056 0.0031 0.0115 0.9998 6493
Colombia 0.0014 0.0003 0.0059 1.0006 5267
Costa Rica 0.0007 0.0001 0.0039 1.0001 4550
Croatia 0.0007 0.0005 0.0011 1.0079 746
Cuba 0.0005 0.0001 0.0035 1.0006 3681
Cyprus 0.0009 0.0007 0.0016 1.0014 1477
Denmark 0.0020 0.0018 0.0024 1.0001 7777
Djibouti 0.0016 0.0009 0.0051 1.0005 3562
Dominican Republic 0.0020 0.0011 0.0057 1.0007 3770
Ecuador 0.0036 0.0019 0.0084 1.0000 7030
Equatorial Guinea 0.0022 0.0006 0.0071 0.9997 8705
Estonia 0.0015 0.0012 0.0020 1.0015 4070
Finland 0.0016 0.0012 0.0031 1.0011 2640
France 0.0023 0.0020 0.0026 0.9996 8243
French Polynesia 0.0008 0.0002 0.0042 1.0000 5227
Gabon 0.0015 0.0004 0.0056 1.0000 7797
Georgia 0.0007 0.0001 0.0045 0.9998 5996
Germany 0.0030 0.0023 0.0051 0.9999 3259
Greece 0.0009 0.0003 0.0044 1.0003 4596
Guinea 0.0012 0.0003 0.0054 1.0000 6866
Hungary 0.0009 0.0004 0.0041 1.0005 3510
Iceland 0.0053 0.0047 0.0060 1.0001 10955
Iran 0.0015 0.0013 0.0020 1.0016 3296
Ireland 0.0050 0.0046 0.0056 1.0000 12822
Israel 0.0023 0.0020 0.0028 1.0002 3992
Italy 0.0042 0.0040 0.0046 0.9997 8790
Jamaica 0.0009 0.0002 0.0046 1.0000 6739
Japan 0.0006 0.0001 0.0037 1.0002 4512
Kazakhstan 0.0009 0.0003 0.0045 0.9997 3775
Kuwait 0.0054 0.0030 0.0112 0.9999 6482
Kyrgyzstan 0.0009 0.0001 0.0048 1.0001 6508
Latvia 0.0013 0.0006 0.0048 0.9998 4046
Lebanon 0.0003 0.0000 0.0031 1.0016 3611
Lithuania 0.0011 0.0006 0.0043 1.0011 4152
Luxembourg 0.0076 0.0061 0.0103 1.0001 4342
Malaysia 0.0010 0.0003 0.0047 1.0003 4832
Maldives 0.0034 0.0014 0.0083 1.0001 6662
Malta 0.0015 0.0010 0.0032 1.0007 2502
Mauritius 0.0007 0.0002 0.0041 1.0003 4226
Mexico 0.0012 0.0004 0.0055 1.0007 4658
Moldova 0.0022 0.0015 0.0043 1.0003 2840
Montenegro 0.0007 0.0004 0.0027 1.0010 2236
Morocco 0.0007 0.0002 0.0045 1.0000 5152
Netherlands 0.0028 0.0026 0.0030 1.0007 10099
New Zealand 0.0012 0.0003 0.0049 0.9999 5150
North Macedonia 0.0008 0.0006 0.0016 1.0022 1610
Norway 0.0020 0.0015 0.0037 1.0015 3039
Oman 0.0015 0.0007 0.0048 0.9998 3873
Palestine 0.0008 0.0001 0.0046 0.9996 6757
Panama 0.0034 0.0025 0.0059 1.0005 3385
Peru 0.0040 0.0022 0.0090 1.0000 5638
Poland 0.0011 0.0006 0.0027 0.9998 2105
Portugal 0.0032 0.0027 0.0045 1.0017 3751
Puerto Rico 0.0017 0.0009 0.0049 1.0015 3759
Qatar 0.0092 0.0066 0.0141 1.0004 7465
Romania 0.0015 0.0010 0.0030 1.0007 2674
Russia 0.0039 0.0023 0.0086 1.0003 4849
Saudi Arabia 0.0014 0.0012 0.0022 1.0005 2893
Serbia 0.0013 0.0012 0.0015 0.9997 5665
Singapore 0.0043 0.0029 0.0076 1.0000 6000
Slovakia 0.0005 0.0002 0.0037 1.0000 4054
Slovenia 0.0008 0.0007 0.0011 1.0008 2591
Spain 0.0052 0.0048 0.0057 0.9997 13208
Sweden 0.0029 0.0026 0.0033 1.0004 8292
Switzerland 0.0039 0.0034 0.0048 0.9998 8101
Turkey 0.0018 0.0018 0.0020 0.9997 10242
Ukraine 0.0013 0.0004 0.0051 1.0001 6063
United Arab Emirates 0.0034 0.0021 0.0074 0.9997 4939
United Kingdom 0.0041 0.0037 0.0048 0.9997 5192
United States 0.0052 0.0047 0.0061 1.0001 5205
Uruguay 0.0004 0.0002 0.0029 1.0001 2498
Beta Albania 0.0512 0.0222 0.0958 1.0004 2729
Armenia 0.0438 0.0261 0.0687 1.0002 7754
Australia 0.0966 0.0838 0.1102 0.9997 10434
Austria 0.1052 0.0895 0.1199 0.9998 6747
Azerbaijan 0.0331 0.0148 0.0652 1.0003 2965
Bahamas 0.0509 0.0216 0.1109 0.9998 4389
Bahrain 0.0689 0.0494 0.0916 0.9998 9855
Barbados 0.0932 0.0514 0.1458 0.9999 8043
Belarus 0.0849 0.0444 0.1459 0.9999 7487
Belgium 0.1166 0.1026 0.1292 0.9996 12681
Bosnia And Herzegovina 0.0421 0.0197 0.0765 0.9999 3785
Brazil 0.0835 0.0616 0.1179 0.9998 6194
Bulgaria 0.0591 0.0222 0.1520 1.0002 3274
Canada 0.0576 0.0392 0.0737 0.9998 6331
Cape Verde 0.0577 0.0259 0.1193 1.0006 7126
Chile 0.0458 0.0361 0.0592 0.9997 10133
Colombia 0.0737 0.0353 0.1249 1.0000 7244
Costa Rica 0.0577 0.0330 0.0881 1.0001 7780
Croatia 0.1051 0.0794 0.1266 1.0010 3542
Cuba 0.0570 0.0243 0.1223 1.0005 4273
Cyprus 0.1086 0.0692 0.1376 1.0005 3111
Denmark 0.0793 0.0639 0.0930 1.0002 7835
Djibouti 0.1527 0.0629 0.2616 1.0007 3608
Dominican Republic 0.0517 0.0278 0.0887 0.9999 6077
Ecuador 0.0823 0.0459 0.1322 1.0001 8194
Equatorial Guinea 0.0757 0.0321 0.1743 0.9997 8599
Estonia 0.0985 0.0710 0.1221 1.0007 5460
Finland 0.0587 0.0283 0.0887 1.0007 3803
France 0.1054 0.0874 0.1219 0.9998 7567
French Polynesia 0.0737 0.0370 0.1241 1.0000 7147
Gabon 0.0695 0.0279 0.1604 0.9999 7444
Georgia 0.0564 0.0231 0.1293 1.0002 6097
Germany 0.0809 0.0563 0.1051 1.0004 4880
Greece 0.0687 0.0393 0.1102 1.0000 5386
Guinea 0.0658 0.0260 0.1546 0.9999 7558
Hungary 0.0530 0.0295 0.0949 0.9998 4087
Iceland 0.1770 0.1626 0.1914 1.0001 10417
Iran 0.0701 0.0465 0.0894 1.0009 4734
Ireland 0.1139 0.0999 0.1271 1.0002 13118
Israel 0.1119 0.0874 0.1339 1.0001 5264
Italy 0.0783 0.0699 0.0862 0.9999 9784
Jamaica 0.0767 0.0269 0.2070 1.0002 5653
Japan 0.0761 0.0323 0.1361 1.0002 6528
Kazakhstan 0.0510 0.0218 0.1108 1.0000 4563
Kuwait 0.0652 0.0488 0.0854 1.0001 9558
Kyrgyzstan 0.0571 0.0239 0.1149 1.0001 7056
Latvia 0.0511 0.0305 0.0831 1.0003 4208
Lebanon 0.0709 0.0255 0.1425 1.0014 5076
Lithuania 0.0751 0.0404 0.1257 1.0007 4917
Luxembourg 0.1067 0.0847 0.1283 1.0000 6311
Malaysia 0.0464 0.0312 0.0683 1.0004 4857
Maldives 0.0930 0.0434 0.2096 1.0002 6137
Malta 0.0862 0.0486 0.1206 1.0011 4075
Mauritius 0.2032 0.1389 0.2903 1.0006 5940
Mexico 0.0550 0.0274 0.1137 1.0004 4789
Moldova 0.0580 0.0286 0.0898 1.0006 4087
Montenegro 0.1500 0.0966 0.2021 1.0009 3586
Morocco 0.0487 0.0205 0.1022 1.0006 5702
Netherlands 0.0982 0.0867 0.1088 1.0001 9709
New Zealand 0.1314 0.1080 0.1580 0.9999 9216
North Macedonia 0.1112 0.0557 0.1537 1.0015 2609
Norway 0.0764 0.0512 0.1002 1.0005 4608
Oman 0.0602 0.0265 0.1243 0.9999 5049
Palestine 0.0623 0.0239 0.1452 1.0000 7850
Panama 0.0554 0.0328 0.0780 1.0000 4909
Peru 0.0781 0.0448 0.1273 0.9998 8437
Poland 0.0372 0.0176 0.0645 0.9999 2825
Portugal 0.0924 0.0618 0.1201 1.0014 4497
Puerto Rico 0.0483 0.0229 0.0923 0.9999 5232
Qatar 0.1023 0.0854 0.1200 0.9997 12195
Romania 0.0500 0.0238 0.0756 1.0006 3659
Russia 0.0824 0.0645 0.1067 0.9999 6254
Saudi Arabia 0.1215 0.0932 0.1435 0.9998 5410
Serbia 0.1136 0.0910 0.1331 0.9997 6483
Singapore 0.1276 0.0870 0.1686 1.0000 6863
Slovakia 0.1255 0.0648 0.2156 1.0001 3519
Slovenia 0.0991 0.0768 0.1183 1.0001 4717
Spain 0.1156 0.1033 0.1267 0.9997 13101
Sweden 0.0893 0.0760 0.1016 1.0003 8691
Switzerland 0.1102 0.0909 0.1277 0.9999 7789
Turkey 0.1138 0.1040 0.1231 0.9999 8660
Ukraine 0.0611 0.0265 0.1309 1.0000 6877
United Arab Emirates 0.0649 0.0351 0.1095 0.9998 6683
United Kingdom 0.0815 0.0663 0.0948 0.9996 6158
United States 0.0773 0.0628 0.0909 0.9997 6312
Uruguay 0.0412 0.0180 0.0819 0.9998 2871
Delta_t Albania 85.8039 45.8600 197.6558 1.0005 2421
Armenia 100.2392 63.9599 168.2192 1.0002 7614
Australia 45.4781 39.8866 52.4659 0.9997 10070
Austria 41.7922 36.6406 49.1236 0.9999 6170
Azerbaijan 132.8678 67.4365 297.2276 1.0004 2536
Bahamas 86.3413 39.6248 203.0833 1.0005 3960
Bahrain 63.7703 47.9797 88.9060 1.0001 8661
Barbados 47.1351 30.1501 85.4934 1.0002 6141
Belarus 51.7464 30.1132 99.0743 0.9999 5947
Belgium 37.6729 34.0154 42.8114 0.9996 11606
Bosnia And Herzegovina 104.3325 57.4332 222.8118 1.0001 3337
Brazil 52.6054 37.2629 71.3522 0.9999 5985
Bulgaria 74.2933 28.9192 197.6372 1.0000 3437
Canada 76.2874 59.5879 112.1083 1.0002 4649
Cape Verde 76.2212 36.8240 169.7384 1.0003 6502
Chile 95.9188 74.1832 121.6169 0.9997 10044
Colombia 59.6663 35.1896 124.3883 1.0002 5213
Costa Rica 76.2139 49.8861 133.3016 1.0006 5980
Croatia 41.7970 34.7211 55.3313 1.0020 2395
Cuba 77.0861 35.9222 181.1727 1.0009 3679
Cyprus 40.4821 31.9312 63.5340 1.0009 2109
Denmark 55.4245 47.2344 68.7461 1.0004 6657
Djibouti 28.7800 16.8011 69.8247 1.0014 3443
Dominican Republic 85.0101 49.5196 158.0127 1.0004 5532
Ecuador 53.3935 33.2344 95.8092 1.0001 6339
Equatorial Guinea 58.0824 25.2149 137.0678 1.0000 9405
Estonia 44.6313 35.9782 61.8825 1.0010 4357
Finland 74.8463 49.5231 155.4443 1.0011 2720
France 41.6769 36.0604 50.2883 0.9998 6670
French Polynesia 59.6650 35.4173 118.9247 1.0001 5978
Gabon 63.1866 27.4052 157.7447 1.0004 6611
Georgia 77.9018 33.9863 190.1147 0.9997 5596
Germany 54.3296 41.8263 78.0767 1.0002 4076
Greece 63.9873 39.8813 111.8177 1.0002 5257
Guinea 66.8117 28.4175 168.9759 0.9999 6291
Hungary 82.8976 46.3219 149.1138 0.9999 3855
Iceland 24.8277 22.9586 27.0182 1.0001 10378
Iran 62.7199 49.1474 94.5981 1.0014 3469
Ireland 38.5897 34.5692 43.9916 1.0003 11561
Israel 39.2661 32.8258 50.2607 1.0001 4189
Italy 56.1221 50.9775 62.8359 0.9999 9276
Jamaica 57.3241 21.2328 163.1457 1.0007 4419
Japan 57.7628 32.2988 136.2348 1.0005 4179
Kazakhstan 86.2233 39.6596 201.2541 1.0002 4282
Kuwait 67.4287 51.4855 90.0861 1.0000 8996
Kyrgyzstan 77.0257 38.2519 183.7266 1.0003 5661
Latvia 85.9917 52.8714 144.0451 1.0002 4089
Lebanon 61.9580 30.8487 172.3209 1.0032 2327
Lithuania 58.4926 34.9689 108.6971 1.0012 4292
Luxembourg 41.1833 34.2462 51.8974 1.0001 5593
Malaysia 94.7399 64.3115 140.8481 1.0005 5128
Maldives 47.2361 20.9646 101.3306 1.0003 6867
Malta 50.9915 36.4511 90.4366 1.0013 3157
Mauritius 21.6276 15.1399 31.6360 1.0003 5104
Mexico 79.9514 38.6394 160.6024 1.0002 4602
Moldova 75.8109 48.9585 153.5525 1.0007 3290
Montenegro 29.2921 21.7491 45.4928 1.0010 2933
Morocco 90.2636 43.0085 214.2565 1.0003 4389
Netherlands 44.7631 40.4014 50.7038 1.0002 9020
New Zealand 33.4484 27.8143 40.6901 1.0000 8904
North Macedonia 39.5124 28.5940 78.9558 1.0021 1829
Norway 57.4830 43.8383 85.7974 1.0010 3799
Oman 72.9524 35.3656 165.9052 0.9997 4714
Palestine 70.4869 30.2626 183.8115 0.9999 6578
Panama 79.2927 56.3702 133.9110 1.0002 3796
Peru 56.2829 34.5199 98.0429 0.9998 7081
Poland 118.0086 68.0864 249.9073 0.9999 2275
Portugal 47.5538 36.5865 71.1074 1.0019 3656
Puerto Rico 90.9290 47.6216 191.6778 1.0005 4795
Qatar 42.9718 36.6290 51.4490 0.9997 11389
Romania 87.8833 58.1402 184.3060 1.0006 2811
Russia 53.3404 41.1924 68.1327 1.0000 5873
Saudi Arabia 36.1545 30.6174 47.1372 0.9999 4073
Serbia 38.6952 33.0223 48.2978 0.9997 5541
Singapore 34.4430 26.0650 50.5321 1.0002 5412
Slovakia 35.0201 20.3778 67.8059 1.0001 3646
Slovenia 44.3220 37.1480 57.2136 1.0002 3712
Spain 38.0207 34.6730 42.5245 0.9996 12597
Sweden 49.2294 43.2520 57.7933 1.0004 7863
Switzerland 39.8783 34.4052 48.3304 1.0000 7087
Turkey 38.6120 35.6971 42.2619 1.0000 8329
Ukraine 71.9490 33.5783 165.7094 0.9997 6562
United Arab Emirates 67.6916 40.1275 125.0253 1.0000 5992
United Kingdom 53.9182 46.3735 66.2336 0.9997 5192
United States 56.8537 48.3290 70.0022 0.9999 5416
Uruguay 106.6129 53.6809 244.2577 1.0004 2671
Alpha Albania -0.1751 -1.2418 1.8334 1.0003 2499
Armenia -2.2600 -3.2804 -1.2989 1.0001 7896
Australia 2.3572 1.0773 3.5762 0.9999 6069
Austria -0.7056 -1.2806 0.0474 1.0000 5769
Azerbaijan -0.5793 -2.4959 1.1442 1.0002 3814
Bahamas -0.0163 -2.8690 2.6180 1.0017 2985
Bahrain -3.6332 -4.3361 -2.9103 0.9998 11971
Barbados 1.4889 -0.4221 3.1425 0.9999 4506
Belarus -2.1618 -3.2868 -0.8582 1.0002 8595
Belgium -2.8410 -3.1796 -2.4524 0.9997 12586
Bosnia And Herzegovina -0.7144 -1.7745 0.7212 1.0003 6021
Brazil -3.0133 -4.0359 -2.4140 1.0006 5262
Bulgaria -1.0753 -3.1048 1.1565 1.0000 4407
Canada -1.8073 -2.0746 -1.5283 0.9996 9482
Cape Verde -1.0277 -3.2384 1.5001 1.0008 4660
Chile -2.8914 -3.6021 -2.3269 0.9998 8053
Colombia -2.8966 -4.5105 -1.0366 1.0010 4686
Costa Rica 2.3286 0.3169 4.0744 1.0004 3627
Croatia -0.6186 -1.2185 0.3885 1.0023 2043
Cuba 1.1499 -2.1020 3.4837 1.0020 2653
Cyprus -1.2195 -1.9779 0.1266 1.0008 2381
Denmark -2.4048 -2.8229 -1.9553 1.0002 8667
Djibouti -0.8632 -2.7144 1.1290 1.0005 3740
Dominican Republic -1.8797 -2.8553 -0.9678 0.9999 8337
Ecuador -3.1784 -4.3966 -1.7246 1.0000 7428
Equatorial Guinea -0.8200 -3.3200 1.9601 1.0000 7364
Estonia -1.5905 -2.3083 -0.6127 1.0009 4871
Finland -1.3304 -2.0343 -0.3633 1.0005 4558
France -2.0570 -2.5405 -1.4683 0.9999 7317
French Polynesia 1.5697 -0.6748 3.4210 1.0001 3774
Gabon -1.2999 -3.7954 2.4193 1.0008 3482
Georgia -2.0003 -4.7451 2.7744 1.0028 2759
Germany -0.7135 -1.5068 0.3569 1.0002 4063
Greece 1.2945 -0.6606 3.1458 1.0001 3457
Guinea -1.5683 -4.1598 2.7903 1.0000 2920
Hungary 0.5736 -0.6426 2.3171 1.0003 3575
Iceland -3.8431 -4.2294 -3.4494 1.0000 10181
Iran -1.7001 -2.3381 -0.8593 1.0007 4640
Ireland -3.0784 -3.3990 -2.7269 1.0002 13277
Israel -1.6297 -2.2173 -0.8765 1.0004 4539
Italy -1.8864 -2.1570 -1.5879 1.0000 9838
Jamaica 1.2548 -3.4010 3.6303 1.0029 2104
Japan 2.4306 -0.2417 4.4559 1.0012 2429
Kazakhstan -1.1201 -3.2699 1.1888 0.9997 5083
Kuwait -2.9839 -3.7297 -2.3161 1.0001 9024
Kyrgyzstan -2.8623 -5.0044 1.9023 1.0009 2568
Latvia 0.5844 -0.8165 2.2099 1.0000 3776
Lebanon 2.7288 -3.0232 5.3177 1.0047 1041
Lithuania 0.0954 -1.5013 1.9343 1.0008 4183
Luxembourg -1.2897 -2.0071 -0.4657 1.0002 5227
Malaysia 1.8346 0.2177 3.4877 1.0001 3792
Maldives -1.3957 -3.1674 0.6982 1.0003 6322
Malta -1.1271 -2.1536 0.4629 1.0007 3313
Mauritius 0.9759 -1.2342 3.0188 1.0013 3818
Mexico -1.8044 -3.5088 -0.3824 1.0002 5095
Moldova -1.1269 -1.6958 -0.3883 1.0003 5984
Montenegro -0.7106 -1.8842 1.5306 1.0011 2523
Morocco -1.8213 -4.2484 2.3382 1.0011 3970
Netherlands -2.0443 -2.3310 -1.7198 1.0001 9833
New Zealand 2.1400 0.5575 3.6223 1.0000 4369
North Macedonia -1.6935 -2.6199 -0.1015 1.0017 2254
Norway -0.9677 -1.8591 0.3243 1.0008 3858
Oman -0.7664 -2.4977 1.1693 1.0001 5309
Palestine -2.9882 -5.5890 3.6356 1.0017 1410
Panama -1.7279 -2.1578 -1.2092 1.0000 10257
Peru -2.3679 -3.3613 -1.2162 1.0001 9169
Poland -0.2013 -0.8306 0.8794 0.9997 3215
Portugal -1.5521 -2.2711 -0.6576 1.0012 4437
Puerto Rico -1.2556 -2.7424 0.4552 1.0001 6443
Qatar -5.1703 -5.7139 -4.5684 0.9999 13994
Romania -0.9905 -1.5026 -0.3074 1.0005 5071
Russia -2.7546 -3.4011 -2.3723 1.0007 5873
Saudi Arabia -2.6758 -2.9019 -2.4389 1.0002 7986
Serbia -1.8155 -2.2058 -1.3528 0.9997 6690
Singapore -4.2295 -5.0474 -3.2884 1.0006 8228
Slovakia 0.4113 -1.6569 2.9171 1.0006 2793
Slovenia -1.2764 -1.8614 -0.4490 1.0002 4033
Spain -2.3982 -2.7346 -2.0275 0.9997 14086
Sweden -3.1338 -3.4458 -2.8092 0.9998 10875
Switzerland -1.7802 -2.3428 -1.1181 0.9998 6475
Turkey -1.9395 -2.1381 -1.7209 1.0003 8954
Ukraine -1.1906 -3.7998 1.9750 1.0003 3686
United Arab Emirates -1.8918 -2.8151 -0.7988 1.0000 7740
United Kingdom -2.1630 -2.4428 -1.8492 0.9996 6979
United States -2.0715 -2.3799 -1.7501 0.9999 7033
Uruguay 0.2195 -1.1643 2.8015 0.9999 2592
Phi Albania 0.1654 0.0809 0.3242 1.0000 12570
Armenia 0.4116 0.2550 0.6811 0.9999 13968
Australia 0.2171 0.1339 0.3660 0.9999 12985
Austria 0.1560 0.1055 0.2430 0.9998 9596
Azerbaijan 0.0737 0.0334 0.1650 1.0003 11287
Bahamas 0.3014 0.0027 1.0863 0.9999 6847
Bahrain 0.6174 0.4224 0.9257 0.9999 13208
Barbados 0.8558 0.2322 2.2635 1.0002 11103
Belarus 0.8221 0.4921 1.4478 0.9996 11999
Belgium 0.0991 0.0669 0.1523 0.9999 11797
Bosnia And Herzegovina 0.3654 0.2235 0.6231 0.9999 12034
Brazil 0.0571 0.0333 0.1109 1.0001 9339
Bulgaria 0.2322 0.1190 0.4797 1.0000 9129
Canada 0.0363 0.0235 0.0604 0.9998 11388
Cape Verde 1.3299 0.6877 2.7268 1.0001 12535
Chile 0.0828 0.0533 0.1390 1.0002 10987
Colombia 0.0699 0.0295 0.2099 0.9999 6586
Costa Rica 0.2913 0.1220 0.6302 0.9997 12376
Croatia 0.0978 0.0512 0.1822 1.0001 9269
Cuba 0.1237 0.0480 0.3314 0.9996 7532
Cyprus 0.2084 0.1092 0.3844 0.9997 8574
Denmark 0.1388 0.0942 0.2121 0.9996 12106
Djibouti 1.1320 0.6674 1.9753 1.0000 8219
Dominican Republic 0.3190 0.1963 0.5515 0.9997 12156
Ecuador 1.8224 1.2460 2.7851 0.9998 11829
Equatorial Guinea 7.2171 3.4994 16.4929 1.0002 8268
Estonia 0.3511 0.2288 0.5547 0.9998 11482
Finland 0.2643 0.1686 0.4272 1.0002 12267
France 0.1809 0.1245 0.2738 0.9998 10606
French Polynesia 0.9928 0.1367 2.9448 0.9998 9470
Gabon 1.5860 0.5663 4.6187 1.0003 8463
Georgia 0.5862 0.2653 1.3374 1.0003 9279
Germany 0.2916 0.1929 0.4620 0.9999 11468
Greece 1.0013 0.6293 1.6477 0.9998 13065
Guinea 1.3362 0.5059 3.9289 1.0008 6791
Hungary 0.1055 0.0600 0.2028 0.9997 11218
Iceland 0.1146 0.0678 0.1999 0.9999 13486
Iran 0.2838 0.1949 0.4277 0.9995 10514
Ireland 0.0965 0.0648 0.1522 0.9998 10417
Israel 0.1992 0.1326 0.3190 0.9996 10520
Italy 0.0520 0.0370 0.0768 0.9998 11958
Jamaica 0.5382 0.1798 1.4986 1.0011 5677
Japan 0.1166 0.0514 0.3371 1.0001 4206
Kazakhstan 0.1465 0.0712 0.3468 0.9997 7911
Kuwait 0.1507 0.0906 0.2672 0.9999 10621
Kyrgyzstan 0.3498 0.1156 1.0112 0.9999 8970
Latvia 0.4341 0.2615 0.7300 0.9997 12784
Lebanon 0.0161 0.0001 0.2232 1.0015 4951
Lithuania 0.8185 0.5180 1.3567 0.9999 12597
Luxembourg 0.3717 0.2498 0.5693 0.9999 12540
Malaysia 0.1514 0.0902 0.2673 0.9996 11826
Maldives 1.3181 0.6841 2.6924 0.9999 11054
Malta 0.5585 0.3382 0.9311 0.9998 11138
Mauritius 0.9286 0.4041 2.1072 0.9998 9435
Mexico 0.0333 0.0146 0.0951 0.9997 7189
Moldova 0.1033 0.0621 0.1826 0.9997 11376
Montenegro 0.3767 0.1695 0.7947 0.9998 8223
Morocco 0.1413 0.0612 0.3840 0.9998 6471
Netherlands 0.0541 0.0360 0.0853 0.9997 11390
New Zealand 0.2856 0.1416 0.5641 1.0000 12172
North Macedonia 0.4509 0.2692 0.7804 1.0002 8913
Norway 0.4475 0.2973 0.6974 1.0000 11995
Oman 0.7307 0.4167 1.3442 0.9997 11616
Palestine 1.5113 0.3787 5.5950 1.0000 7793
Panama 0.1065 0.0694 0.1715 0.9996 10726
Peru 0.5795 0.3540 1.0218 0.9997 11047
Poland 0.0411 0.0239 0.0748 0.9999 11652
Portugal 0.2849 0.1887 0.4518 0.9998 12302
Puerto Rico 0.9210 0.5744 1.5606 0.9997 11440
Qatar 0.5272 0.3576 0.7939 0.9998 13978
Romania 0.0743 0.0472 0.1237 0.9997 11686
Russia 0.0212 0.0121 0.0426 0.9997 8375
Saudi Arabia 0.0225 0.0127 0.0454 1.0000 7780
Serbia 0.0614 0.0372 0.1090 0.9999 10018
Singapore 0.8508 0.5541 1.3544 0.9999 12348
Slovakia 0.6435 0.3704 1.1348 0.9998 8405
Slovenia 0.1714 0.0911 0.3197 0.9998 9311
Spain 0.0910 0.0621 0.1397 0.9998 11065
Sweden 0.0857 0.0576 0.1314 0.9998 11815
Switzerland 0.2539 0.1722 0.3895 1.0002 11378
Turkey 0.0188 0.0119 0.0316 1.0002 11395
Ukraine 1.8238 1.0519 3.3788 1.0000 11824
United Arab Emirates 0.5374 0.3290 0.9212 0.9996 11510
United Kingdom 0.0445 0.0296 0.0718 1.0000 10085
United States 0.0502 0.0334 0.0798 0.9997 11604
Uruguay 0.1949 0.0869 0.3954 0.9998 14560 
p <- results %>% 
  ungroup %>% 
  filter(par %in% c("beta", "alpha", "S", "phi_inv", "nu", "delta_t")) %>% 
  mutate(par = str_replace(par, "phi_inv", "phi")) %>% 
  inner_join(countries, by = c("num" = "country_id")) %>% 
  mutate(par = str_to_title(par)) %>% 
  mutate(plot_var = str_c(par, "_", country)) %>% 
  ggplot(aes(country, estimate, ymin = conf.low, ymax = conf.high, col = country == "Iceland")) +
  geom_linerange() +
  geom_point() +
  coord_flip() +
  facet_wrap("par", scales = "free_x", ncol = 2) +
  scale_colour_manual(values = c("grey", "blue")) +
  labs(title = "Posterior medians and 95% PIs") +
  theme(axis.title = element_blank())

ggplotly(p)

Hyperparameters

results %>% 
  ungroup %>% 
  filter(par %in% c("mu_beta", "sigma_beta", 
                    "mu_delta_t",
                    "mu_alpha", "sigma_alpha", 
                    "mu_s", "kappa_s",
                    "sigma_phi_inv_sqrt", "mu_nu", "sigma_nu")) %>% 
  mutate(par = str_replace(par, "sigma_phi_inv_sqrt", "sigma_phi_sqrt")) %>% 
  select(par, everything(), -num, -term, -std.error) %>% 
  set_names(c("Parameter",  "Median", "Lower", "Upper", "Rhat", "ESS")) %>% 
  kable(digits = 4, align = c("l", "l", rep("c", ncol(.) - 2)),
        caption = "Table 3. Summary of posterior samples of hyperparameters") %>% 
  kable_styling(bootstrap_options = c("striped", "hover")) %>% 
  add_header_above(c("", "", "95% PI" = 2, "Convergence" = 2)) %>% 
  column_spec(1, bold = T)
Table 3. Summary of posterior samples of hyperparameters
95% PI
Convergence
Parameter Median Lower Upper Rhat ESS
mu_beta -2.5793 -2.7383 -2.4519 1.0032 973
sigma_beta 0.4571 0.3474 0.6055 1.0005 1202
mu_alpha -1.1826 -1.6272 -0.7136 1.0054 629
sigma_alpha 1.8184 1.4463 2.2892 1.0013 1175
mu_s 0.0025 0.0020 0.0031 1.0001 3137
kappa_s 637.0741 415.3739 967.3564 1.0003 3079
sigma_phi_sqrt 0.6069 0.4903 0.7635 1.0077 319
mu_delta_t 57.9540 51.0188 67.9395 1.0033 977

Posterior Predictions

We obtain predictions by conditioning on the first observed number of cumulative cases in a country (as seen in the modeling data), perform simulations of new daily cases and sum those up to get predicted cumulative cases.

Posterior predictions for \(R_e\) are obtained by looking at part of the SIR differential equation on daily new cases

\[ \frac{d}{dt} I = \beta \cdot I \cdot \frac{S}{N} - \text{Some rate of recovery} \]

We do not model the rate of recovery, but based on literature review we assume that diagnosed cases are healthy after 15 days. As we are not interested in the negative part of \(\frac{d}{dt}I\) we write \(\frac{d}{dt} I_{i, t} = D_{i, t}\), the number of new cases. We can then simplify to estimate \(\beta\) in country \(i\) at time \(t\) as

\[ \beta_{i, t} = D_{i, t} \cdot \frac{N}{I_{i, t} \cdot S_{i, t}} \]

where \(I_{i, t}\) is the cumulative number of diagnosed cases at time \(t\) minus the cumulative number of cases at time \(t - 15\). \(S_{i, t}\) is simply the population minus the cumulative number of cases at time \(t\). Then we utilize

\[ R_e = \beta \gamma^{-1} \]

where \(\gamma^{-1}\) is the length of the serial interval. After a literature review we put a \(\mathrm{Gamma}(30, 2.5)\) distibution on the length of serial interval.

knitr::include_app("https://bgautijonsson.shinyapps.io/Hierarchical_Report_NegBin/", height = "900px")

Validation

GIFs with historical predictions for some countries made by the model can be seen here.

Some things to keep in mind

Comparisons between countries

This model does not predict the total amount of COVID-19 cases, it predicts the total amount of DIAGNOSED COVID-19 cases. It then depends on each country’s specific diagnostic criteria what that number means. Iceland does many tests per capita so there it is predicting a broader spectrum of cases, not all serious, whereas in a country like Italy or Spain the model’s predictions have another meaning as the diagnostic criteria are different from those of Iceland.

Changes in policies

As of yet, the model only implements one slope for each country, which is assumed to be fixed throughout the epidemic. Each country’s policies should of course affect this growth rate and the date at which a country implements a policy should affect when the slope changes. Finding solutions to this is on our next actions list.

Next actions

  • Time-varying growth rates.
  • Include country-level information in the model

Estimating the burden on the Icelandic health care system

We use the current age distribution of cases in Iceland to split our posterior predictive distritubion of daily cases into age categories. We then use numbers from Ferguson et al. at Imperial College, (Table 1), to sample numbers of patients in hospital and from those numbers, how many will end up in ICU. Along with the numbers above, we use the following parameters in our simulations

  • lag_infected_to_hospital <- 7
  • lag_hospital_to_icu <- 3
  • days_from_infection_to_healthy <- 21
  • days_in_hospital <- 14
  • days_in_icu <- 10

The results from our predictions are then shared at covid.hi.is

Appendix

Stan Code

data {
  int<lower = 0> N_obs;
  // country[N_obs] takes on value i for country number i etc.
  int country[N_obs];
  // days is how many days since case rate per 1000 exceeded a limit we have chosen
  vector[N_obs] days;
  // We model the daily number of newly diagnosed cases instead of cumulative numbers
  int new_cases[N_obs];
  
  int<lower = 0> N_countries;
  vector[N_countries] pop;
}

parameters {
  // Since we use a non-centered parametrisation we first create normal(0, 1) variables for alpha and beta
  
  // Beta parameters
  vector[N_countries] z_beta;
  real mu_beta;
  real<lower = 0> sigma_beta;
  
  // Alpha parameters
  vector[N_countries] z_alpha;
  real mu_alpha;
  real<lower = 0> sigma_alpha;
  
  //  Asymptote/saturation parameters
  // We model the beta distribution in terms of mean and sample size instead of a and b.
  vector<lower = 0, upper = 1>[N_countries] S;
  real<lower = 0, upper = 1> mu_s;
  real<lower = 0> kappa_s;
  
  //  Overdispersion parameters. One for each country which is dependent on hyperparameter
  vector<lower = 0>[N_countries] z_phi_inv_sqrt;
  real<lower = 0> sigma_phi_inv_sqrt;
  
}

transformed parameters {
   // Non-Centerd parametrisations
   // If B ~ normal(mu_b, sigma_b) then B = mu_b + sigma_b * normal(0, 1)
  vector<lower = 0>[N_countries] beta = exp(mu_beta + sigma_beta * z_beta);
  vector[N_countries] alpha = mu_alpha + sigma_alpha * z_alpha;
   // If X ~ exponential(lambda) then X ~ lambda * exponential(1)
  vector<lower = 0>[N_countries] phi_inv_sqrt = sigma_phi_inv_sqrt * z_phi_inv_sqrt;
   // Overdispersion parameters
  vector<lower = 0>[N_countries] phi_inv = square(phi_inv_sqrt);
  vector<lower = 0>[N_countries] phi = inv(phi_inv);
   // Asymptote hyperparameters
  real<lower = 0> a_s = mu_s * kappa_s;
  real<lower = 0> b_s = (1 - mu_s) * kappa_s;
   // Logistic equation calculations
  vector[N_obs] linear = alpha[country] + beta[country] .* days;
  vector<lower = 0>[N_obs] difference;
  for (i in 1:N_obs) {
    difference[i] = beta[country[i]] * S[country[i]] * exp(-linear[i]) / square(exp(-linear[i]) + 1);
  }
}

model {
   // Alpha parameters
  mu_alpha ~ normal(-2.5, 3);
  sigma_alpha ~ exponential(1);
  z_alpha ~ std_normal();
  
   // Beta parameters
  mu_beta ~ normal(-3, 1);
  sigma_beta ~ exponential(1);
  z_beta ~ std_normal();
  
   // Asymptote parameters
  mu_s ~ beta(1, 99);
  kappa_s ~ exponential(0.001);
  S ~ beta(a_s, b_s);
  
   // Overdispersion parameters
  z_phi_inv_sqrt ~ exponential(1);
  sigma_phi_inv_sqrt ~ exponential(1);
  
  //  Likelihood
  new_cases ~ neg_binomial_2(difference .* pop[country], phi[country]);
}