assess new `socialmixr` interface after 0.6.0 CRAN release

[avallecam]

Currently opting to keep workflow 2, unless demography gets embedded in workflow 3

#pak::pak("epiverse-trace/epidemics@fix/model-default-demography-index")

# 1 old --------------------
# data(polymod, package = "socialmixr")

# contact_data1 <- socialmixr::contact_matrix(
#   polymod,
#   countries = "United Kingdom",
#   age_limits = c(0, 20, 40),
#   symmetric = TRUE,
#   return_demography = TRUE
# )

# contact_matrix1 <- t(contact_data1$matrix)
# demography_vector1 <- contact_data1$demography$population
# names(demography_vector1) <- rownames(contact_matrix1)


# 2 new: using survey_pop ------------------------------------------------------------
data(polymod, package = "socialmixr")
data(popAge1dt, package = "wpp2024")

uk_pop <- popAge1dt |>
  dplyr::filter(name == "United Kingdom", year == max(year)) |>
  dplyr::select(lower.age.limit = age, population = pop) |>
  dplyr::mutate(population = population * 1000)

contact_data2 <- socialmixr::contact_matrix(
  polymod,
  countries = "United Kingdom",
  age_limits = c(0, 20, 40),
  symmetric = TRUE,
  survey_pop = uk_pop
)

contact_matrix2 <- t(contact_data2$matrix)
demography_vector2 <- contact_data2$demography$population
names(demography_vector2) <- rownames(contact_matrix2)


# 3 new: using tidyverse ------------------------------------------------------------
data(polymod, package = "socialmixr")
data(popAge1dt, package = "wpp2024")

uk_pop <- popAge1dt |> 
  dplyr::filter(name == "United Kingdom", year == max(year)) |> 
  dplyr::select(lower.age.limit = age, population = pop) |> 
  dplyr::mutate(population = population * 1000)

# given that we need it twice
uk_age_limits <- c(0, 20, 40)

contact_data3 <- polymod[country == "United Kingdom"] |>
  socialmixr::assign_age_groups(age_limits = uk_age_limits) |>
  socialmixr::compute_matrix() |>
  socialmixr::symmetrise(survey_pop = uk_pop) # request

contact_matrix3 <- t(contact_data3$matrix)

# given that there is no 
# return_demography = TRUE
demography_vector3 <- uk_pop |>
  dplyr::mutate(
    age_category = base::cut(
      x = lower.age.limit,
      breaks = c(uk_age_limits, Inf), # requires explicit Inf
      include.lowest = TRUE,
      right = FALSE
    )
  ) |> 
  dplyr::group_by(age_category) |> 
  dplyr::summarise(population = sum(population)) |> 
  dplyr::pull(population, name = "age_category")


# outputs ------------------------------------------------------------

#contact_matrix1
contact_matrix2
#>                  age.group
#> contact.age.group   [0,20)  [20,40) [40,Inf)
#>          [0,20)   7.883663 2.806513 1.488256
#>          [20,40)  3.105519 4.854839 2.495156
#>          [40,Inf) 3.290997 4.986325 5.005571
contact_matrix3
#>                  age.group
#> contact.age.group   [0,20)  [20,40) [40,Inf)
#>          [0,20)   7.883663 2.806513 1.488256
#>          [20,40)  3.105519 4.854839 2.495156
#>          [40,Inf) 3.290997 4.986325 5.005571

#demography_vector1
demography_vector2
#>   [0,20)  [20,40) [40,Inf) 
#> 15962345 17662976 35297722
demography_vector3
#>   [0,20)  [20,40) [40,Inf] 
#> 15962345 17662976 35297722

contact_matrix <- contact_matrix3
demography_vector <- demography_vector3


# downstream --------------------------------------------------------

# initial conditions
initial_i <- 1e-6
initial_conditions <- c(
  S = 1 - initial_i, E = 0, I = initial_i, R = 0, V = 0
)

# build for all age groups
initial_conditions <- rbind(
  initial_conditions,
  initial_conditions,
  initial_conditions
)

# assign rownames for clarity
rownames(initial_conditions) <- rownames(contact_matrix)

# view initial conditions
initial_conditions
#>                 S E     I R V
#> [0,20)   0.999999 0 1e-06 0 0
#> [20,40)  0.999999 0 1e-06 0 0
#> [40,Inf) 0.999999 0 1e-06 0 0

# population object
uk_population <- epidemics::population(
  name = "UK",
  contact_matrix = contact_matrix,
  demography_vector = demography_vector,
  initial_conditions = initial_conditions
)

uk_population
#> <population> object
#> 
#>  Population name:
#> "UK"
#> 
#>  Demography 
#> [0,20): 15,962,345 (20%)
#> [20,40): 17,662,976 (30%)
#> [40,Inf): 35,297,722 (50%)
#> 
#>  Contact matrix 
#>                  age.group
#> contact.age.group   [0,20)  [20,40) [40,Inf)
#>          [0,20)   7.883663 2.806513 1.488256
#>          [20,40)  3.105519 4.854839 2.495156
#>          [40,Inf) 3.290997 4.986325 5.005571
#> 
#>  Initial Conditions 
#>                 S E     I R V
#> [0,20)   0.999999 0 1e-06 0 0
#> [20,40)  0.999999 0 1e-06 0 0
#> [40,Inf) 0.999999 0 1e-06 0 0

^{Created on 2026-06-05 with reprex v2.1.1}

[avallecam]

note:

• assesss if max(year) is appropriate to keep