Final size attack rate simulation from contact data

finalsize_contact.R

# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
# Illustrative code to general finalsize model outputs from contact matrix
# Adapted from: https://epiverse-trace.github.io/finalsize/articles/finalsize.html
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 

library(devtools)
install.packages("finalsize")
library(finalsize)
library(rio)
library(tidyverse)
remotes::install_cran("socialmixr"); library(socialmixr)

age_labels <- seq(0,80,10)

# Load UK social mixing data ----------------------------------------------

contact_data <- contact_matrix(
  polymod,
  countries = "United Kingdom",
  age.limits = age_labels, # Use age bands
  symmetric = TRUE
)

# Define model parameters -------------------------------------------------

demography_vector = contact_data$demography$population
contact_matrix = t(contact_data$matrix)
contact_matrix = contact_matrix / max(eigen(contact_matrix)$values)

names(demography_vector) = contact_data$demography$age.group
n_demo_grps = length(demography_vector)

contact_matrix = contact_matrix / demography_vector

# Define susceptibility based on proportion infected
susceptibility = matrix(
  data = 1, # Susceptibility same across groups (i.e. new virus)
  n_demo_grps
)

n_risk_grps = 1L
p_susceptibility = matrix(1, n_demo_grps, n_risk_grps)


# Run model
r0 = 2

output_model <- finalsize::final_size(
  r0 = r0,
  contact_matrix = contact_matrix,
  demography_vector = demography_vector,
  susceptibility = susceptibility,
  p_susceptibility = p_susceptibility,
  solver = "iterative"
)


plot(tail(age_labels,-1),output_model$p_infected,col="blue",pch=19,lwd=2,
     xlab="age",
     ylab="proportion infected",
     ylim=c(0,1))