Transmission Metrics

Calculate transmission metrics to estimate epidemic potential, velocity, and spread

Effective Reproduction

The Effective Reproduction Number (Re) is the average number of secondary premises infected by a source per day. Re is used to estimate epidemic potential.

Code
## exclude index cases (source_farm == 0)
transmissions <- merge %>% 
  filter(source_farm != 0)

## Count the number of transmissions per source farm for each iteration on each infect_day
iteration_Re <- transmissions %>%
  group_by(region, scenario_type, preclinical, iteration, infect_day, source_farm) %>%
  summarize(num_transmissions = n(), .groups = "drop") 

## Filter to regions
iteration_Re_western <- iteration_Re %>%
  filter(region == "western")

iteration_Re_central <- iteration_Re %>%
  filter(region == "central")

iteration_Re_eastern <- iteration_Re %>%
  filter(region == "eastern")

View example of daily iteration_Re

Code
## Filter to scenario to display
iteration_Re_central_select <- iteration_Re_central %>%
  filter(scenario_type == "suboptimal") %>%
  filter(preclinical == "2")

## Select and order columns to display
iteration_Re_central_select <- iteration_Re_central_select[c("iteration", "infect_day", "source_farm", "num_transmissions")]

head(iteration_Re_central_select)
iteration infect_day source_farm num_transmissions
1 14 788420 3
1 15 788420 4
1 16 788088 1
1 16 788420 1
1 17 788088 2
1 19 788420 1

Daily Summaries

Calculate daily summaries for each region.

Code
## Western U.S.
daily_Re_western <- calculate_daily_Re(iteration_Re_western)

## Central U.S.
daily_Re_central <- calculate_daily_Re(iteration_Re_central)

## Eastern U.S.
daily_Re_eastern <- calculate_daily_Re(iteration_Re_eastern)

Significance Testing

Perform significance testing on optimal and suboptimal detection scenarios.

Linear Model

Code
no_LV_western_Re <-daily_Re_western %>%
  filter(scenario_type != "low-virulence")

model_western_Re <- lm(q50 ~ preclinical * scenario_type, data = no_LV_western_Re)

summary(model_western_Re)

Call:
lm(formula = q50 ~ preclinical * scenario_type, data = no_LV_western_Re)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.0279 -0.3610 -0.0279  0.3551  2.6854 

Coefficients:
                             Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   0.41468    0.01782  23.272  < 2e-16 ***
preclinical1                  0.08052    0.02520   3.195  0.00141 ** 
preclinical2                  0.10218    0.02520   4.055 5.15e-05 ***
preclinical3                  0.31592    0.02520  12.537  < 2e-16 ***
scenario_type.L               0.14152    0.02520   5.616 2.14e-08 ***
preclinical1:scenario_type.L  0.07018    0.03564   1.969  0.04902 *  
preclinical2:scenario_type.L  0.07895    0.03564   2.215  0.02681 *  
preclinical3:scenario_type.L  0.27894    0.03564   7.827 6.99e-15 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.4755 on 2840 degrees of freedom
Multiple R-squared:  0.1813,    Adjusted R-squared:  0.1793 
F-statistic: 89.85 on 7 and 2840 DF,  p-value: < 2.2e-16
Code
no_LV_central_Re <-daily_Re_central %>%
  filter(scenario_type != "low-virulence")

model_central_Re <- lm(q50 ~ preclinical * scenario_type, data = no_LV_central_Re)

summary(model_central_Re)

Call:
lm(formula = q50 ~ preclinical * scenario_type, data = no_LV_central_Re)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.0248 -0.4965 -0.0248  0.3988  3.6667 

Coefficients:
                             Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   0.41491    0.01805  22.987  < 2e-16 ***
preclinical1                  0.08392    0.02553   3.288  0.00102 ** 
preclinical2                  0.21631    0.02553   8.474  < 2e-16 ***
preclinical3                  0.39809    0.02553  15.595  < 2e-16 ***
scenario_type.L               0.11537    0.02553   4.520 6.45e-06 ***
preclinical1:scenario_type.L  0.15642    0.03610   4.333 1.52e-05 ***
preclinical2:scenario_type.L  0.05658    0.03610   1.567  0.11713    
preclinical3:scenario_type.L  0.18420    0.03610   5.102 3.57e-07 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.4816 on 2840 degrees of freedom
Multiple R-squared:  0.1731,    Adjusted R-squared:  0.171 
F-statistic: 84.91 on 7 and 2840 DF,  p-value: < 2.2e-16
Code
no_LV_eastern_Re <-daily_Re_eastern %>%
  filter(scenario_type != "low-virulence")

model_eastern_Re <- lm(q50 ~ preclinical * scenario_type, data = no_LV_eastern_Re)

summary(model_eastern_Re)

Call:
lm(formula = q50 ~ preclinical * scenario_type, data = no_LV_eastern_Re)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.97114 -0.35112  0.02886  0.35896  2.75772 

Coefficients:
                             Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   0.40403    0.01732  23.323  < 2e-16 ***
preclinical1                  0.13970    0.02450   5.702 1.30e-08 ***
preclinical2                  0.38666    0.02450  15.783  < 2e-16 ***
preclinical3                  0.29502    0.02450  12.042  < 2e-16 ***
scenario_type.L               0.22876    0.02450   9.337  < 2e-16 ***
preclinical1:scenario_type.L  0.04363    0.03465   1.259    0.208    
preclinical2:scenario_type.L -0.01713    0.03465  -0.494    0.621    
preclinical3:scenario_type.L  0.15603    0.03465   4.504 6.95e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.4622 on 2840 degrees of freedom
Multiple R-squared:  0.2248,    Adjusted R-squared:  0.2228 
F-statistic: 117.6 on 7 and 2840 DF,  p-value: < 2.2e-16

ANOVA

Code
anova(model_western_Re)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 3 38.85586 12.9519547 57.29261 0
scenario_type 1 87.96228 87.9622845 389.09873 0
preclinical:scenario_type 3 15.36251 5.1208378 22.65188 0
Residuals 2840 642.02956 0.2260667 NA NA
Code
anova(model_central_Re)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 3 64.360300 21.4534335 92.48313 0e+00
scenario_type 1 65.621438 65.6214383 282.88601 0e+00
preclinical:scenario_type 3 7.887333 2.6291108 11.33378 2e-07
Residuals 2840 658.798525 0.2319713 NA NA
Code
anova(model_eastern_Re)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 3 62.222792 20.7409306 97.07196 0.0e+00
scenario_type 1 107.212903 107.2129031 501.77915 0.0e+00
preclinical:scenario_type 3 6.483673 2.1612242 10.11499 1.3e-06
Residuals 2840 606.810079 0.2136655 NA NA

Plot Effective Reproduction

The dashed line is at Re = 1

Code
central_Re_plot <- plot_daily_Re(daily_Re_central)

central_Re_plot

Code
eastern_Re_plot <- plot_daily_Re(daily_Re_eastern)

eastern_Re_plot

Epidemic Velocity

The epidemic velocity is the distance (km) of spread per day. The epidemic velocity is estimated by first calculating the daily spread distances.

Calculate Daily Distances

compile_daily_summary() calculates distances (km) between source farms and those infected for each iteration, then calculates the percentiles by day across all iterations to get the average statistics.

Code
daily_summary_western <- compile_daily_summary(western_paths, reference)

daily_summary_central <- compile_daily_summary(central_paths, reference)

daily_summary_eastern <- compile_daily_summary(eastern_paths, reference)

Plot Epidemic Velocity

Code
## Get quantiles to compare scenarios
daily_central <- daily_summary_central$combined_summary

## Plot epidemic velocity
central_velocity_plot <- plot_wave_velocity(daily_central)

central_velocity_plot

Code
## Get quantiles to compare scenarios
daily_eastern <- daily_summary_eastern$combined_summary

## Plot epidemic velocity
eastern_velocity_plot <- plot_wave_velocity(daily_eastern)

eastern_velocity_plot

Iteration Metrics

Compare scenarios by looking for patterns at the level of individual iterations. iteration_metrics() will use the daily_distances output to look at trends between optimal and suboptimal responses.

The function returns iteration specific metrics:

auc_log

The amount of area (total) under the plotted line (distance x day) above on a log scale, i.e., the Area Under Curve (AUC). This represents the total distance covered during the outbreak`s spread (on the log scale).

peak_spread

The maximum distance spread in any one day (log scale).

peak_day

The infect day that the maximum spread distance (peak_spread) occurred.

Code
## Calculate daily iteration metrics
iteration_metrics_western <- iteration_metrics(daily_summary_western$daily_distances) %>%
  select(-c(scenario))

iteration_metrics_central <- iteration_metrics(daily_summary_central$daily_distances) %>%
  select(-c(scenario))

iteration_metrics_eastern <- iteration_metrics(daily_summary_eastern$daily_distances) %>%
  select(-c(scenario))

View example of iteration metrics

Code
## Filter data to scenario
iteration_metrics_central_select <- iteration_metrics_central %>%
  filter(scenario_type == "Suboptimal Detection") %>%
  filter(preclinical == "3")

## Select and order columns to display
iteration_metrics_central_select <- iteration_metrics_central_select[c("iteration", "auc_log", "peak_spread", "peak_day")]

head(iteration_metrics_central_select)
iteration auc_log peak_spread peak_day
1 770.6477 6.177126 64
2 502.7660 5.571266 49
3 100.8983 4.881051 17
4 355.1653 6.656522 62
5 726.0179 6.111867 18
6 182.3753 5.375928 43

Cumulative Spread

Drop the low-virulence scenarios to compare optimal and suboptimal detection scenarios. auc_log represents the total cumulative spread distance, for each scenario by iteration, on a log scale.

Code
## Drop the low-virulence scenario
iteration_metrics_no_low_western <- iteration_metrics_western %>%
  filter(scenario_type != "Low-Virulence")

## Plot AUC
ggplot(iteration_metrics_no_low_western, 
       aes(x = preclinical, y = auc_log, color = scenario_type)) +
  geom_point(shape=1, alpha = 0.3, size = 2) + # individual iterations
  geom_smooth(method = "lm", se = FALSE, linewidth=1.2) + # trend
  ylim(0, 1500) +
  facet_grid(. ~ scenario_type) +
  scale_color_manual(values = c("Suboptimal Detection" = "#74add1", "Optimal Detection" = "orange2")) +
  labs(
    x = "Preclinical Infectious Duration (days)",
    y = "Cumulative Spread (log AUC)",
    title = " "
  ) +
  theme_minimal() +
  theme(
    plot.margin    = unit(c(0.25, 0.25, 0.25, 0.25), "cm"),
    legend.position = "none",
    strip.text     = element_text(size = 18, face = "bold", color = "gray40"),
    axis.title.x   = element_text(size = 22, face = "bold"),
    axis.title.y   = element_text(size = 22, face = "bold"),
    axis.text.x    = element_text(size = 18, face = "bold"),
    axis.text.y    = element_text(size = 18, face = "bold"),
    plot.title     = element_text(size = 22, face = "bold", hjust = 0.5)
  )

Code
## Drop the low-virulence scenario
iteration_metrics_no_low_central <- iteration_metrics_central %>%
  filter(scenario_type != "Low-Virulence")

## Plot AUC
ggplot(iteration_metrics_no_low_central, 
       aes(x = preclinical, y = auc_log, color = scenario_type)) +
  geom_point(shape=1, alpha = 0.3, size = 2) + # individual iterations
  geom_smooth(method = "lm", se = FALSE, linewidth=1.2) + # trend
  ylim(0, 1500) +
  facet_grid(. ~ scenario_type) +
  scale_color_manual(values = c("Suboptimal Detection" = "#74add1", "Optimal Detection" = "orange2")) +
  labs(
    x = "Preclinical Infectious Duration (days)",
    y = "Cumulative Spread (log AUC)",
    title = " "
  ) +
  theme_minimal() +
  theme(
    plot.margin    = unit(c(0.25, 0.25, 0.25, 0.25), "cm"),
    legend.position = "none",
    strip.text     = element_text(size = 18, face = "bold", color = "gray40"),
    axis.title.x   = element_text(size = 22, face = "bold"),
    axis.title.y   = element_text(size = 22, face = "bold"),
    axis.text.x    = element_text(size = 18, face = "bold"),
    axis.text.y    = element_text(size = 18, face = "bold"),
    plot.title     = element_text(size = 22, face = "bold", hjust = 0.5)
  )

Code
## Drop the low-virulence scenario
iteration_metrics_no_low_eastern <- iteration_metrics_eastern %>%
  filter(scenario_type != "Low-Virulence")

## Plot AUC
ggplot(iteration_metrics_no_low_eastern, 
       aes(x = preclinical, y = auc_log, color = scenario_type)) +
  geom_point(shape=1, alpha = 0.3, size = 2) + # individual iterations
  geom_smooth(method = "lm", se = FALSE, linewidth=1.2) + # trend
  ylim(0, 1500) +
  facet_grid(. ~ scenario_type) +
  scale_color_manual(values = c("Suboptimal Detection" = "#74add1", "Optimal Detection" = "orange2")) +
  labs(
    x = "Preclinical Infectious Duration (days)",
    y = "Cumulative Spread (log AUC)",
    title = " "
  ) +
  theme_minimal() +
  theme(
    plot.margin    = unit(c(0.25, 0.25, 0.25, 0.25), "cm"),
    legend.position = "none",
    strip.text     = element_text(size = 18, face = "bold", color = "gray40"),
    axis.title.x   = element_text(size = 22, face = "bold"),
    axis.title.y   = element_text(size = 22, face = "bold"),
    axis.text.x    = element_text(size = 18, face = "bold"),
    axis.text.y    = element_text(size = 18, face = "bold"),
    plot.title     = element_text(size = 22, face = "bold", hjust = 0.5)
  )

Significance Test

Fit a linear model with an interaction term to quantify the influence of preclinical infectious duration (preclinical) and detection scenario (scenario_type) on cumulative outbreak spread (auc_log).

Linear Model

Code
model_auc_western <- lm(auc_log ~ preclinical * scenario_type, data = iteration_metrics_no_low_western)

summary(model_auc_western)

Call:
lm(formula = auc_log ~ preclinical * scenario_type, data = iteration_metrics_no_low_western)

Residuals:
    Min      1Q  Median      3Q     Max 
-283.79  -73.51  -21.68   59.97  957.56 

Coefficients:
                            Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   94.283      3.569  26.415   <2e-16 ***
preclinical                   38.029      1.824  20.846   <2e-16 ***
scenario_type.L               54.450      5.048  10.787   <2e-16 ***
preclinical:scenario_type.L   22.173      2.580   8.594   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 117.4 on 3500 degrees of freedom
  (277 observations deleted due to missingness)
Multiple R-squared:  0.2976,    Adjusted R-squared:  0.297 
F-statistic: 494.2 on 3 and 3500 DF,  p-value: < 2.2e-16
Code
model_auc_central <- lm(auc_log ~ preclinical * scenario_type, data = iteration_metrics_no_low_central)

summary(model_auc_central)

Call:
lm(formula = auc_log ~ preclinical * scenario_type, data = iteration_metrics_no_low_central)

Residuals:
    Min      1Q  Median      3Q     Max 
-379.62  -92.16  -32.04   64.68  918.16 

Coefficients:
                            Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   79.564      5.467  14.554  < 2e-16 ***
preclinical                   58.007      2.753  21.070  < 2e-16 ***
scenario_type.L               55.472      7.731   7.175 8.85e-13 ***
preclinical:scenario_type.L   42.508      3.893  10.918  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 172.8 on 3391 degrees of freedom
  (301 observations deleted due to missingness)
Multiple R-squared:  0.294, Adjusted R-squared:  0.2934 
F-statistic: 470.8 on 3 and 3391 DF,  p-value: < 2.2e-16
Code
model_auc_eastern <- lm(auc_log ~ preclinical * scenario_type, data = iteration_metrics_no_low_eastern)

summary(model_auc_eastern)

Call:
lm(formula = auc_log ~ preclinical * scenario_type, data = iteration_metrics_no_low_eastern)

Residuals:
    Min      1Q  Median      3Q     Max 
-316.97  -79.18  -28.81   47.52  855.07 

Coefficients:
                            Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   78.708      4.674  16.840  < 2e-16 ***
preclinical                   45.739      2.373  19.276  < 2e-16 ***
scenario_type.L               48.654      6.610   7.361 2.28e-13 ***
preclinical:scenario_type.L   32.891      3.356   9.802  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 150.3 on 3429 degrees of freedom
  (280 observations deleted due to missingness)
Multiple R-squared:  0.2637,    Adjusted R-squared:  0.263 
F-statistic: 409.3 on 3 and 3429 DF,  p-value: < 2.2e-16

ANOVA

Code
anova(model_auc_western)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 1 5445104 5445103.59 395.05040 0
scenario_type 1 13972417 13972417.05 1013.71972 0
preclinical:scenario_type 1 1018076 1018076.45 73.86297 0
Residuals 3500 48241598 13783.31 NA NA
Code
anova(model_auc_central)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 1 12687728 12687727.73 424.8123 0
scenario_type 1 25933670 25933669.92 868.3148 0
preclinical:scenario_type 1 3560019 3560019.23 119.1971 0
Residuals 3391 101277875 29866.67 NA NA
Code
anova(model_auc_eastern)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 1 8014986 8014986.44 354.99578 0
scenario_type 1 17537524 17537523.84 776.76326 0
preclinical:scenario_type 1 2169123 2169123.39 96.07373 0
Residuals 3429 77418915 22577.69 NA NA

Peak Spread

As with cumulative spread above, comparisons can be made using the peak spread rates.

Code
ggplot(iteration_metrics_no_low_western, aes(x = preclinical, y = peak_day, color = scenario_type)) +
  geom_point(shape=1, alpha = 0.3, size = 2) + # individual iterations
  geom_smooth(method = "lm", se = FALSE, linewidth=1.2) +
  ylim(0, 365) +
  facet_grid(. ~ scenario_type) +
  scale_color_manual(values = c("Suboptimal Detection" = "#74add1", "Optimal Detection" = "orange2")) +
  labs(
    x = "Preclinical Infectious Duration (days)",
    y = "Day of Peak Transmission",
    title = " "
  ) +
  theme_minimal() +
  theme(
    plot.margin    = unit(c(0.25, 0.25, 0.25, 0.25), "cm"),
    legend.position = "none",
    strip.text     = element_text(size = 18, face = "bold", color = "gray40"),
    axis.title.x   = element_text(size = 22, face = "bold"),
    axis.title.y   = element_text(size = 22, face = "bold"),
    axis.text.x    = element_text(size = 18, face = "bold"),
    axis.text.y    = element_text(size = 18, face = "bold"),
    plot.title     = element_text(size = 22, face = "bold", hjust = 0.5)
  )

Code
ggplot(iteration_metrics_no_low_central, aes(x = preclinical, y = peak_day, color = scenario_type)) +
  geom_point(shape=1, alpha = 0.3, size = 2) + # individual iterations
  geom_smooth(method = "lm", se = FALSE, linewidth=1.2) +
  ylim(0, 365) +
  facet_grid(. ~ scenario_type) +
  scale_color_manual(values = c("Suboptimal Detection" = "#74add1", "Optimal Detection" = "orange2")) +
  labs(
    x = "Preclinical Infectious Duration (days)",
    y = "Day of Peak Transmission",
    title = " "
  ) +
  theme_minimal() +
  theme(
    plot.margin    = unit(c(0.25, 0.25, 0.25, 0.25), "cm"),
    legend.position = "none",
    strip.text     = element_text(size = 18, face = "bold", color = "gray40"),
    axis.title.x   = element_text(size = 22, face = "bold"),
    axis.title.y   = element_text(size = 22, face = "bold"),
    axis.text.x    = element_text(size = 18, face = "bold"),
    axis.text.y    = element_text(size = 18, face = "bold"),
    plot.title     = element_text(size = 22, face = "bold", hjust = 0.5)
  )

Code
ggplot(iteration_metrics_no_low_eastern, aes(x = preclinical, y = peak_day, color = scenario_type)) +
  geom_point(shape=1, alpha = 0.3, size = 2) + # individual iterations
  geom_smooth(method = "lm", se = FALSE, linewidth=1.2) +
  ylim(0, 365) +
  facet_grid(. ~ scenario_type) +
  scale_color_manual(values = c("Suboptimal Detection" = "#74add1", "Optimal Detection" = "orange2")) +
  labs(
    x = "Preclinical Infectious Duration (days)",
    y = "Day of Peak Transmission",
    title = " "
  ) +
  theme_minimal() +
  theme(
    plot.margin    = unit(c(0.25, 0.25, 0.25, 0.25), "cm"),
    legend.position = "none",
    strip.text     = element_text(size = 18, face = "bold", color = "gray40"),
    axis.title.x   = element_text(size = 22, face = "bold"),
    axis.title.y   = element_text(size = 22, face = "bold"),
    axis.text.x    = element_text(size = 18, face = "bold"),
    axis.text.y    = element_text(size = 18, face = "bold"),
    plot.title     = element_text(size = 22, face = "bold", hjust = 0.5)
  )

Significance Test

Fit a linear model with an interaction term to quantify the influence of preclinical and scenario_type on the maximum spread rate and it’s (peak_day).

Linear Model

Code
model_peak_western <- lm(peak_day ~ preclinical * scenario_type, data = iteration_metrics_no_low_western)

summary(model_peak_western)

Call:
lm(formula = peak_day ~ preclinical * scenario_type, data = iteration_metrics_no_low_western)

Residuals:
    Min      1Q  Median      3Q     Max 
-45.185 -14.874  -6.384   8.801 193.815 

Coefficients:
                            Estimate Std. Error t value Pr(>|t|)    
(Intercept)                  29.3035     0.7215  40.617  < 2e-16 ***
preclinical                   4.6295     0.3783  12.237  < 2e-16 ***
scenario_type.L               8.5991     1.0203   8.428  < 2e-16 ***
preclinical:scenario_type.L   3.2587     0.5350   6.091 1.24e-09 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 25.76 on 3777 degrees of freedom
Multiple R-squared:  0.1571,    Adjusted R-squared:  0.1564 
F-statistic: 234.7 on 3 and 3777 DF,  p-value: < 2.2e-16
Code
model_peak_central <- lm(peak_day ~ preclinical * scenario_type, data = iteration_metrics_no_low_central)

summary(model_peak_central)

Call:
lm(formula = peak_day ~ preclinical * scenario_type, data = iteration_metrics_no_low_central)

Residuals:
    Min      1Q  Median      3Q     Max 
-56.527 -15.082  -6.755   5.828 229.473 

Coefficients:
                            Estimate Std. Error t value Pr(>|t|)    
(Intercept)                  25.1348     1.0158  24.745  < 2e-16 ***
preclinical                   7.5614     0.5279  14.322  < 2e-16 ***
scenario_type.L               7.1520     1.4365   4.979 6.69e-07 ***
preclinical:scenario_type.L   6.9064     0.7466   9.250  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 35.46 on 3692 degrees of freedom
Multiple R-squared:  0.1699,    Adjusted R-squared:  0.1692 
F-statistic: 251.8 on 3 and 3692 DF,  p-value: < 2.2e-16
Code
model_peak_eastern <- lm(peak_day ~ preclinical * scenario_type, data = iteration_metrics_no_low_eastern)

summary(model_peak_eastern)

Call:
lm(formula = peak_day ~ preclinical * scenario_type, data = iteration_metrics_no_low_eastern)

Residuals:
    Min      1Q  Median      3Q     Max 
-43.909 -13.940  -6.919   5.852 287.081 

Coefficients:
                            Estimate Std. Error t value Pr(>|t|)    
(Intercept)                  26.3329     0.8619  30.553  < 2e-16 ***
preclinical                   4.7061     0.4506  10.444  < 2e-16 ***
scenario_type.L               6.5160     1.2189   5.346 9.54e-08 ***
preclinical:scenario_type.L   4.6434     0.6373   7.286 3.87e-13 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 30.38 on 3709 degrees of freedom
Multiple R-squared:  0.1263,    Adjusted R-squared:  0.1256 
F-statistic: 178.7 on 3 and 3709 DF,  p-value: < 2.2e-16

ANOVA

Code
anova(model_peak_western)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 1 91984.87 91984.8676 138.58997 0
scenario_type 1 350638.02 350638.0205 528.29246 0
preclinical:scenario_type 1 24620.58 24620.5842 37.09486 0
Residuals 3777 2506868.63 663.7195 NA NA
Code
anova(model_peak_central)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 1 248069.9 248069.913 197.27732 0
scenario_type 1 594273.2 594273.198 472.59509 0
preclinical:scenario_type 1 107599.6 107599.568 85.56844 0
Residuals 3692 4642571.8 1257.468 NA NA
Code
anova(model_peak_eastern)
Df Sum Sq Mean Sq F value Pr(>F)
preclinical 1 97130.83 97130.8281 105.25678 0
scenario_type 1 348554.42 348554.4173 377.71442 0
preclinical:scenario_type 1 48993.79 48993.7929 53.09261 0
Residuals 3709 3422660.79 922.7988 NA NA