Skip to content

Commit 1c4723f

Browse files
committed
differences for PR #145
1 parent 2e59211 commit 1c4723f

8 files changed

Lines changed: 147 additions & 91 deletions

.DS_Store

-6 KB
Binary file not shown.

compare-interventions.md

Lines changed: 16 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -236,19 +236,25 @@ output <- epidemics::model_vacamole(
236236

237237

238238
``` r
239-
polymod <- socialmixr::polymod
240-
contact_data <- socialmixr::contact_matrix(
241-
survey = polymod,
239+
survey_files_2 <- contactsurveys::download_survey(
240+
survey = "https://doi.org/10.5281/zenodo.3874557",
241+
verbose = FALSE
242+
)
243+
survey_load_2 <- socialmixr::load_survey(files = survey_files_2)
244+
245+
contacts_byage_2 <- socialmixr::contact_matrix(
246+
survey = survey_load_2,
242247
countries = "United Kingdom",
243-
age.limits = c(0, 20, 40),
244-
symmetric = TRUE
248+
age_limits = c(0, 20, 40),
249+
symmetric = TRUE,
250+
return_demography = TRUE
245251
)
246252
# prepare contact matrix
247-
contact_matrix <- t(contact_data$matrix)
253+
contacts_byage_matrix_2 <- t(contacts_byage_2$matrix)
248254

249255
# extract demography vector
250-
demography_vector <- contact_data$demography$population
251-
names(demography_vector) <- rownames(contact_matrix)
256+
demography_vector <- contacts_byage_2$demography$population
257+
names(demography_vector) <- rownames(contacts_byage_matrix_2)
252258

253259
# prepare initial conditions
254260
initial_i <- 1e-6
@@ -266,12 +272,12 @@ initial_conditions_vacamole <- rbind(
266272
initial_conditions_vacamole,
267273
initial_conditions_vacamole
268274
)
269-
rownames(initial_conditions_vacamole) <- rownames(contact_matrix)
275+
rownames(initial_conditions_vacamole) <- rownames(contacts_byage_matrix_2)
270276

271277
# prepare population object
272278
uk_population_vacamole <- epidemics::population(
273279
name = "UK",
274-
contact_matrix = contact_matrix,
280+
contact_matrix = contacts_byage_matrix_2,
275281
demography_vector = demography_vector,
276282
initial_conditions = initial_conditions_vacamole
277283
)

contact-matrices.md

Lines changed: 51 additions & 22 deletions
Original file line numberDiff line numberDiff line change
@@ -36,6 +36,7 @@ Some groups of individuals have more contacts than others; the average schoolchi
3636

3737

3838
``` r
39+
library(contactsurveys)
3940
library(socialmixr)
4041
```
4142

@@ -78,25 +79,49 @@ For a contact matrix with rows $i$ and columns $j$:
7879

7980
Contact matrices are commonly estimated from studies that use diaries to record interactions. For example, the POLYMOD survey measured contact patterns in 8 European countries using data on the location and duration of contacts reported by the study participants [(Mossong et al. 2008)](https://doi.org/10.1371/journal.pmed.0050074).
8081

81-
The R package `{socialmixr}` contains functions which can estimate contact matrices from POLYMOD and other surveys. We can load the POLYMOD survey data:
82+
The R package `{socialmixr}` contains functions which can estimate contact matrices from POLYMOD and other surveys. We can download and load the POLYMOD survey data directly from Zenodo using `{contactsurveys}` and `{socialmixr}`:
8283

8384

85+
``` r
86+
survey_files <- contactsurveys::download_survey(
87+
survey = "https://doi.org/10.5281/zenodo.3874557",
88+
verbose = FALSE
89+
)
90+
91+
survey_load <- socialmixr::load_survey(files = survey_files)
92+
```
93+
94+
::::::::::::::::::::::::::::::::::::: callout
95+
### Specify the country name
96+
97+
A single survey file can contain data from multiple countries. You can inspect the available countries with:
98+
8499

85100
``` r
86-
polymod <- socialmixr::polymod
101+
levels(survey_load$participants$country)
87102
```
88103

89-
Then we can obtain the contact matrix for the age categories we want by specifying `age_limits`.
104+
``` output
105+
[1] "Belgium" "Finland" "Germany" "Italy"
106+
[5] "Luxembourg" "Netherlands" "Poland" "United Kingdom"
107+
```
108+
109+
Always pass the `countries =` argument to `contact_matrix()` to make sure you use data from the intended country only.
110+
111+
::::::::::::::::::::::::::::::::::::::::::::::::
112+
113+
Then we can obtain the contact matrix for the age categories we want by specifying `age_limits`. We also add `return_demography = TRUE` to include demographic information in the output, which is required when using the contact matrix with `{epidemics}`.
90114

91115

92116
``` r
93-
contact_data <- socialmixr::contact_matrix(
94-
survey = polymod,
117+
contacts_byage <- socialmixr::contact_matrix(
118+
survey = survey_load,
95119
countries = "United Kingdom",
96120
age_limits = c(0, 20, 40),
97-
symmetric = TRUE
121+
symmetric = TRUE,
122+
return_demography = TRUE
98123
)
99-
contact_data
124+
contacts_byage
100125
```
101126

102127
``` output
@@ -124,8 +149,8 @@ $participants
124149

125150

126151

127-
**Note: although the contact matrix `contact_data$matrix` is not itself mathematically symmetric, it satisfies the condition that the total number of contacts of one group with another is the same as the reverse. In other words:
128-
`contact_data$matrix[j,i]*contact_data$demography$proportion[j] = contact_data$matrix[i,j]*contact_data$demography$proportion[i]`.
152+
**Note: although the contact matrix `contacts_byage$matrix` is not itself mathematically symmetric, it satisfies the condition that the total number of contacts of one group with another is the same as the reverse. In other words:
153+
`contacts_byage$matrix[j,i]*contacts_byage$demography$proportion[j] = contacts_byage$matrix[i,j]*contacts_byage$demography$proportion[i]`.
129154
For the mathematical explanation see [the corresponding section in the socialmixr documentation](https://epiforecasts.io/socialmixr/articles/socialmixr.html#symmetric-contact-matrices).**
130155

131156

@@ -145,14 +170,17 @@ If `symmetric` is set to TRUE, the `contact_matrix()` function will internally u
145170

146171
::::::::::::::::::::::::::::::::::::::::::::::::
147172

148-
The example above uses the POLYMOD survey. There are a number of surveys available in `socialmixr`. To list the available surveys, use `socialmixr::list_surveys()`. To download a survey, we can use `socialmixr::get_survey()`
173+
The example above uses the POLYMOD survey. There are a number of surveys available in `socialmixr`. To list the available surveys, use `socialmixr::list_surveys()`. To download a survey from Zenodo and load it, we use `contactsurveys::download_survey()` followed by `socialmixr::load_survey()`:
149174

150175

151176
``` r
152-
# Access the contact survey data from Zenodo
153-
zambia_sa_survey <- socialmixr::get_survey(
154-
"https://doi.org/10.5281/zenodo.3874675"
177+
# Download and load the contact survey data for Zambia from Zenodo
178+
zambia_survey_files <- contactsurveys::download_survey(
179+
survey = "https://doi.org/10.5281/zenodo.3874675",
180+
verbose = FALSE
155181
)
182+
183+
zambia_sa_survey <- socialmixr::load_survey(files = zambia_survey_files)
156184
```
157185

158186
:::::::::::::::::: spoiler
@@ -211,11 +239,12 @@ Similar to the code above, to access vector values within a dataframe, you can u
211239

212240
``` r
213241
# Generate the contact matrix for Zambia only
214-
contact_data_zambia <- socialmixr::contact_matrix(
242+
contacts_byage_zambia <- socialmixr::contact_matrix(
215243
survey = zambia_sa_survey,
216244
countries = "Zambia", # key argument
217245
age_limits = c(0, 20),
218-
symmetric = TRUE
246+
symmetric = TRUE,
247+
return_demography = TRUE
219248
)
220249
```
221250

@@ -228,7 +257,7 @@ participants).
228257

229258
``` r
230259
# Print the contact matrix for Zambia only
231-
contact_data_zambia
260+
contacts_byage_zambia
232261
```
233262

234263
``` output
@@ -253,7 +282,7 @@ $participants
253282

254283
``` r
255284
# Print the vector of population size for {epidemics}
256-
contact_data_zambia$demography$population
285+
contacts_byage_zambia$demography$population
257286
```
258287

259288
``` output
@@ -329,13 +358,13 @@ When simulating an epidemic, we often want to ensure that the average number of
329358

330359
Rather than just using the raw number of contacts, we can instead normalise the contact matrix to make it easier to work in terms of $R_0$. In particular, we normalise the matrix by scaling it so that if we were to calculate the average number of secondary cases based on this normalised matrix, the result would be 1 (in mathematical terms, we are scaling the matrix so the largest eigenvalue is 1). This transformation scales the entries but preserves their relative values.
331360

332-
In the case of the above model, we want to define $\beta C_{i,j}$ so that the model has a specified valued of $R_0$. If the entry of the contact matrix $C[i,j]$ represents the contacts of population $i$ with $j$, it is equivalent to `contact_data$matrix[i,j]`, and the maximum eigenvalue of this matrix represents the typical magnitude of contacts, not typical magnitude of transmission. We must therefore normalise the matrix $C$ so the maximum eigenvalue is one; we call this matrix $C_{normalised}$. Because the rate of recovery is $\gamma$, individuals will be infectious on average for $1/\gamma$ days. So $\beta$ as a model input is calculated from $R_0$, the scaling factor and the value of $\gamma$ (i.e. mathematically we use the fact that the dominant eigenvalue of the matrix $R_0 \times C_{normalised}$ is equal to $\beta / \gamma$).
361+
In the case of the above model, we want to define $\beta C_{i,j}$ so that the model has a specified valued of $R_0$. If the entry of the contact matrix $C[i,j]$ represents the contacts of population $i$ with $j$, it is equivalent to `contacts_byage$matrix[i,j]`, and the maximum eigenvalue of this matrix represents the typical magnitude of contacts, not typical magnitude of transmission. We must therefore normalise the matrix $C$ so the maximum eigenvalue is one; we call this matrix $C_{normalised}$. Because the rate of recovery is $\gamma$, individuals will be infectious on average for $1/\gamma$ days. So $\beta$ as a model input is calculated from $R_0$, the scaling factor and the value of $\gamma$ (i.e. mathematically we use the fact that the dominant eigenvalue of the matrix $R_0 \times C_{normalised}$ is equal to $\beta / \gamma$).
333362

334363

335364
``` r
336-
contact_matrix <- t(contact_data$matrix)
337-
scaling_factor <- 1 / max(eigen(contact_matrix)$values)
338-
normalised_matrix <- contact_matrix * scaling_factor
365+
contacts_byage_matrix <- t(contacts_byage$matrix)
366+
scaling_factor <- 1 / max(eigen(contacts_byage_matrix)$values)
367+
normalised_matrix <- contacts_byage_matrix * scaling_factor
339368
```
340369

341370
As a result, if we multiply the scaled matrix by $R_0$, then converting to the number of expected secondary cases would give us $R_0$, as required.
@@ -363,7 +392,7 @@ Normalisation can be performed by the function `contact_matrix()` in `{socialmix
363392

364393
``` r
365394
contact_data_split <- socialmixr::contact_matrix(
366-
survey = polymod,
395+
survey = survey_load,
367396
countries = "United Kingdom",
368397
age_limits = c(0, 20, 40),
369398
symmetric = TRUE,

disease-burden.md

Lines changed: 17 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -48,6 +48,7 @@ In this tutorial, we'll focus on the separated approach, where we first run an e
4848
``` r
4949
library(epiparameter)
5050
library(epidemics)
51+
library(contactsurveys)
5152
library(socialmixr)
5253
library(tidyverse)
5354
```
@@ -85,21 +86,27 @@ We'll use `{epiparameter}` to define these delay distributions. The Gamma distri
8586

8687

8788
``` r
88-
# load contact and population data from socialmixr::polymod
89-
polymod <- socialmixr::polymod
90-
contact_data <- socialmixr::contact_matrix(
91-
polymod,
89+
# download and load contact and population data from Zenodo
90+
survey_files <- contactsurveys::download_survey(
91+
survey = "https://doi.org/10.5281/zenodo.3874557",
92+
verbose = FALSE
93+
)
94+
survey_load <- socialmixr::load_survey(files = survey_files)
95+
96+
contacts_byage <- socialmixr::contact_matrix(
97+
survey = survey_load,
9298
countries = "United Kingdom",
9399
age_limits = c(0, 20, 40),
94-
symmetric = TRUE
100+
symmetric = TRUE,
101+
return_demography = TRUE
95102
)
96103

97104
# prepare contact matrix
98-
contact_matrix <- t(contact_data$matrix)
105+
contacts_byage_matrix <- t(contacts_byage$matrix)
99106

100107
# prepare the demography vector
101-
demography_vector <- contact_data$demography$population
102-
names(demography_vector) <- rownames(contact_matrix)
108+
demography_vector <- contacts_byage$demography$population
109+
names(demography_vector) <- rownames(contacts_byage_matrix)
103110

104111
# initial conditions: one in every 1 million is infected
105112
initial_i <- 1e-6
@@ -117,12 +124,12 @@ initial_conditions <- rbind(
117124
initial_conditions_free,
118125
initial_conditions_free
119126
)
120-
rownames(initial_conditions) <- rownames(contact_matrix)
127+
rownames(initial_conditions) <- rownames(contacts_byage_matrix)
121128

122129
# prepare the population to model as affected by the epidemic
123130
uk_population <- epidemics::population(
124131
name = "UK",
125-
contact_matrix = contact_matrix,
132+
contact_matrix = contacts_byage_matrix,
126133
demography_vector = demography_vector,
127134
initial_conditions = initial_conditions
128135
)

md5sum.txt

Lines changed: 7 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -4,16 +4,16 @@
44
"config.yaml" "a0c04c1d43ce0640c3ea333b140e89c8" "site/built/config.yaml" "2026-03-20"
55
"index.md" "32bc80d6f4816435cc0e01540cb2a513" "site/built/index.md" "2026-03-20"
66
"links.md" "8184cf4149eafbf03ce8da8ff0778c14" "site/built/links.md" "2026-03-20"
7-
"episodes/contact-matrices.Rmd" "747723602382243f24cf566ea2d344d0" "site/built/contact-matrices.md" "2026-03-20"
8-
"episodes/simulating-transmission.Rmd" "9d5d9a20894160f25d93bddb761bd2ee" "site/built/simulating-transmission.md" "2026-03-20"
7+
"episodes/contact-matrices.Rmd" "0bcf84fd375c77d20b0391baebe826fb" "site/built/contact-matrices.md" "2026-03-20"
8+
"episodes/simulating-transmission.Rmd" "ac558c97836c06a076f28d2d6a469765" "site/built/simulating-transmission.md" "2026-03-20"
99
"episodes/model-choices.Rmd" "aa195e66455fb6a97b4930fd08c08001" "site/built/model-choices.md" "2026-03-20"
10-
"episodes/modelling-interventions.Rmd" "2cce81c1ceb8ad0caef9b9e0e584e753" "site/built/modelling-interventions.md" "2026-03-20"
11-
"episodes/compare-interventions.Rmd" "5675995a87e94e383a5a138a57e0f906" "site/built/compare-interventions.md" "2026-03-20"
12-
"episodes/vaccine-comparisons.Rmd" "b59e2eed4984431a20d1f9b3d8f3fe58" "site/built/vaccine-comparisons.md" "2026-03-20"
13-
"episodes/disease-burden.Rmd" "501ccbed5022e40c25ee467cd97fead1" "site/built/disease-burden.md" "2026-03-20"
10+
"episodes/modelling-interventions.Rmd" "e9501b49a94ca9ffa515cd95c94c6872" "site/built/modelling-interventions.md" "2026-03-20"
11+
"episodes/compare-interventions.Rmd" "3f995f8fcdd60fa3b57ddbcbe7b6bc9f" "site/built/compare-interventions.md" "2026-03-20"
12+
"episodes/vaccine-comparisons.Rmd" "a31706c4185b80c29fe05f5549b55ca5" "site/built/vaccine-comparisons.md" "2026-03-20"
13+
"episodes/disease-burden.Rmd" "65768cdde110fa85019d8e31256790b1" "site/built/disease-burden.md" "2026-03-20"
1414
"instructors/instructor-notes.md" "ca3834a1b0f9e70c4702aa7a367a6bb5" "site/built/instructor-notes.md" "2026-03-20"
1515
"learners/BF_measles.Rmd" "5725a25d664730bbce3736a3243c2182" "site/built/BF_measles.md" "2026-03-20"
1616
"learners/reference.md" "9e836f1ec999f95f55135cdd2c78d2e6" "site/built/reference.md" "2026-03-20"
1717
"learners/setup.md" "be1de7ec41f52f14d00bcb74e87b5046" "site/built/setup.md" "2026-03-20"
1818
"profiles/learner-profiles.md" "31b503c4b5bd1f0960ada730eca4a25e" "site/built/learner-profiles.md" "2026-03-20"
19-
"renv/profiles/lesson-requirements/renv.lock" "933ec57112a484dff787ccaf4ac71829" "site/built/renv.lock" "2026-03-20"
19+
"renv/profiles/lesson-requirements/renv.lock" "fc9f27801a86b09d34d57df019048b09" "site/built/renv.lock" "2026-03-20"

modelling-interventions.md

Lines changed: 19 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -55,6 +55,7 @@ In this tutorial, we will learn how to use `{epidemics}` to model interventions
5555

5656
``` r
5757
library(epidemics)
58+
library(contactsurveys)
5859
library(socialmixr)
5960
library(tidyverse)
6061
```
@@ -84,23 +85,28 @@ The SEIR model divides the population into four compartments: Susceptible (S), E
8485

8586

8687
``` r
87-
# load survey data
88-
survey_data <- socialmixr::polymod
88+
# download and load survey data
89+
survey_files <- contactsurveys::download_survey(
90+
survey = "https://doi.org/10.5281/zenodo.3874557",
91+
verbose = FALSE
92+
)
93+
survey_load <- socialmixr::load_survey(files = survey_files)
8994

9095
# generate contact matrix
91-
cm_results <- socialmixr::contact_matrix(
92-
survey = survey_data,
96+
contacts_byage <- socialmixr::contact_matrix(
97+
survey = survey_load,
9398
countries = "United Kingdom",
9499
age_limits = c(0, 15, 65),
95-
symmetric = TRUE
100+
symmetric = TRUE,
101+
return_demography = TRUE
96102
)
97103

98104
# transpose contact matrix
99-
cm_matrix <- t(cm_results$matrix)
105+
contacts_byage_matrix <- t(contacts_byage$matrix)
100106

101107
# prepare the demography vector
102-
demography_vector <- cm_results$demography$population
103-
names(demography_vector) <- rownames(cm_matrix)
108+
demography_vector <- contacts_byage$demography$population
109+
names(demography_vector) <- rownames(contacts_byage_matrix)
104110

105111
# initial conditions: one in every 1 million is infected
106112
initial_i <- 1e-6
@@ -118,12 +124,12 @@ initial_conditions <- base::rbind(
118124
initial_conditions,
119125
initial_conditions
120126
)
121-
rownames(initial_conditions) <- rownames(cm_matrix)
127+
rownames(initial_conditions) <- rownames(contacts_byage_matrix)
122128

123129
# prepare the population to model as affected by the epidemic
124130
uk_population <- epidemics::population(
125131
name = "UK",
126-
contact_matrix = cm_matrix,
132+
contact_matrix = contacts_byage_matrix,
127133
demography_vector = demography_vector,
128134
initial_conditions = initial_conditions
129135
)
@@ -177,7 +183,7 @@ To include an intervention in our model we must create an `intervention` object.
177183

178184

179185
``` r
180-
rownames(cm_matrix)
186+
rownames(contacts_byage_matrix)
181187
```
182188

183189
``` output
@@ -419,8 +425,8 @@ Here we will assume all age groups are vaccinated at the same rate 0.01 and that
419425
# prepare a vaccination object
420426
vaccinate <- epidemics::vaccination(
421427
name = "vaccinate all",
422-
time_begin = matrix(40, nrow(cm_matrix)),
423-
time_end = matrix(40 + 150, nrow(cm_matrix)),
428+
time_begin = matrix(40, nrow(contacts_byage_matrix)),
429+
time_end = matrix(40 + 150, nrow(contacts_byage_matrix)),
424430
nu = matrix(c(0.01, 0.01, 0.01))
425431
)
426432
```

0 commit comments

Comments
 (0)