Skip to content

Commit a8be09b

Browse files
committed
differences for PR #145
1 parent f128305 commit a8be09b

8 files changed

Lines changed: 178 additions & 91 deletions

.DS_Store

-6 KB
Binary file not shown.

compare-interventions.md

Lines changed: 22 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -236,19 +236,31 @@ 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+
246+
``` warning
247+
Warning: No reference provided.
248+
```
249+
250+
``` r
251+
contacts_byage_2 <- socialmixr::contact_matrix(
252+
survey = survey_load_2,
242253
countries = "United Kingdom",
243-
age.limits = c(0, 20, 40),
244-
symmetric = TRUE
254+
age_limits = c(0, 20, 40),
255+
symmetric = TRUE,
256+
return_demography = TRUE
245257
)
246258
# prepare contact matrix
247-
contact_matrix <- t(contact_data$matrix)
259+
contacts_byage_matrix_2 <- t(contacts_byage_2$matrix)
248260

249261
# extract demography vector
250-
demography_vector <- contact_data$demography$population
251-
names(demography_vector) <- rownames(contact_matrix)
262+
demography_vector <- contacts_byage_2$demography$population
263+
names(demography_vector) <- rownames(contacts_byage_matrix_2)
252264

253265
# prepare initial conditions
254266
initial_i <- 1e-6
@@ -266,12 +278,12 @@ initial_conditions_vacamole <- rbind(
266278
initial_conditions_vacamole,
267279
initial_conditions_vacamole
268280
)
269-
rownames(initial_conditions_vacamole) <- rownames(contact_matrix)
281+
rownames(initial_conditions_vacamole) <- rownames(contacts_byage_matrix_2)
270282

271283
# prepare population object
272284
uk_population_vacamole <- epidemics::population(
273285
name = "UK",
274-
contact_matrix = contact_matrix,
286+
contact_matrix = contacts_byage_matrix_2,
275287
demography_vector = demography_vector,
276288
initial_conditions = initial_conditions_vacamole
277289
)

contact-matrices.md

Lines changed: 55 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,53 @@ 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+
``` warning
95+
Warning: No reference provided.
96+
```
97+
98+
::::::::::::::::::::::::::::::::::::: callout
99+
### Specify the country name
100+
101+
A single survey file can contain data from multiple countries. You can inspect the available countries with:
102+
84103

85104
``` r
86-
polymod <- socialmixr::polymod
105+
levels(survey_load$participants$country)
87106
```
88107

89-
Then we can obtain the contact matrix for the age categories we want by specifying `age_limits`.
108+
``` output
109+
[1] "Belgium" "Finland" "Germany" "Italy"
110+
[5] "Luxembourg" "Netherlands" "Poland" "United Kingdom"
111+
```
112+
113+
Always pass the `countries =` argument to `contact_matrix()` to make sure you use data from the intended country only.
114+
115+
::::::::::::::::::::::::::::::::::::::::::::::::
116+
117+
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}`.
90118

91119

92120
``` r
93-
contact_data <- socialmixr::contact_matrix(
94-
survey = polymod,
121+
contacts_byage <- socialmixr::contact_matrix(
122+
survey = survey_load,
95123
countries = "United Kingdom",
96124
age_limits = c(0, 20, 40),
97-
symmetric = TRUE
125+
symmetric = TRUE,
126+
return_demography = TRUE
98127
)
99-
contact_data
128+
contacts_byage
100129
```
101130

102131
``` output
@@ -124,8 +153,8 @@ $participants
124153

125154

126155

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]`.
156+
**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:
157+
`contacts_byage$matrix[j,i]*contacts_byage$demography$proportion[j] = contacts_byage$matrix[i,j]*contacts_byage$demography$proportion[i]`.
129158
For the mathematical explanation see [the corresponding section in the socialmixr documentation](https://epiforecasts.io/socialmixr/articles/socialmixr.html#symmetric-contact-matrices).**
130159

131160

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

146175
::::::::::::::::::::::::::::::::::::::::::::::::
147176

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()`
177+
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()`:
149178

150179

151180
``` r
152-
# Access the contact survey data from Zenodo
153-
zambia_sa_survey <- socialmixr::get_survey(
154-
"https://doi.org/10.5281/zenodo.3874675"
181+
# Download and load the contact survey data for Zambia from Zenodo
182+
zambia_survey_files <- contactsurveys::download_survey(
183+
survey = "https://doi.org/10.5281/zenodo.3874675",
184+
verbose = FALSE
155185
)
186+
187+
zambia_sa_survey <- socialmixr::load_survey(files = zambia_survey_files)
156188
```
157189

158190
:::::::::::::::::: spoiler
@@ -211,11 +243,12 @@ Similar to the code above, to access vector values within a dataframe, you can u
211243

212244
``` r
213245
# Generate the contact matrix for Zambia only
214-
contact_data_zambia <- socialmixr::contact_matrix(
246+
contacts_byage_zambia <- socialmixr::contact_matrix(
215247
survey = zambia_sa_survey,
216248
countries = "Zambia", # key argument
217249
age_limits = c(0, 20),
218-
symmetric = TRUE
250+
symmetric = TRUE,
251+
return_demography = TRUE
219252
)
220253
```
221254

@@ -228,7 +261,7 @@ participants).
228261

229262
``` r
230263
# Print the contact matrix for Zambia only
231-
contact_data_zambia
264+
contacts_byage_zambia
232265
```
233266

234267
``` output
@@ -253,7 +286,7 @@ $participants
253286

254287
``` r
255288
# Print the vector of population size for {epidemics}
256-
contact_data_zambia$demography$population
289+
contacts_byage_zambia$demography$population
257290
```
258291

259292
``` output
@@ -329,13 +362,13 @@ When simulating an epidemic, we often want to ensure that the average number of
329362

330363
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.
331364

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$).
365+
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$).
333366

334367

335368
``` r
336-
contact_matrix <- t(contact_data$matrix)
337-
scaling_factor <- 1 / max(eigen(contact_matrix)$values)
338-
normalised_matrix <- contact_matrix * scaling_factor
369+
contacts_byage_matrix <- t(contacts_byage$matrix)
370+
scaling_factor <- 1 / max(eigen(contacts_byage_matrix)$values)
371+
normalised_matrix <- contacts_byage_matrix * scaling_factor
339372
```
340373

341374
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 +396,7 @@ Normalisation can be performed by the function `contact_matrix()` in `{socialmix
363396

364397
``` r
365398
contact_data_split <- socialmixr::contact_matrix(
366-
survey = polymod,
399+
survey = survey_load,
367400
countries = "United Kingdom",
368401
age_limits = c(0, 20, 40),
369402
symmetric = TRUE,

disease-burden.md

Lines changed: 23 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,33 @@ 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+
97+
``` warning
98+
Warning: No reference provided.
99+
```
100+
101+
``` r
102+
contacts_byage <- socialmixr::contact_matrix(
103+
survey = survey_load,
92104
countries = "United Kingdom",
93105
age_limits = c(0, 20, 40),
94-
symmetric = TRUE
106+
symmetric = TRUE,
107+
return_demography = TRUE
95108
)
96109

97110
# prepare contact matrix
98-
contact_matrix <- t(contact_data$matrix)
111+
contacts_byage_matrix <- t(contacts_byage$matrix)
99112

100113
# prepare the demography vector
101-
demography_vector <- contact_data$demography$population
102-
names(demography_vector) <- rownames(contact_matrix)
114+
demography_vector <- contacts_byage$demography$population
115+
names(demography_vector) <- rownames(contacts_byage_matrix)
103116

104117
# initial conditions: one in every 1 million is infected
105118
initial_i <- 1e-6
@@ -117,12 +130,12 @@ initial_conditions <- rbind(
117130
initial_conditions_free,
118131
initial_conditions_free
119132
)
120-
rownames(initial_conditions) <- rownames(contact_matrix)
133+
rownames(initial_conditions) <- rownames(contacts_byage_matrix)
121134

122135
# prepare the population to model as affected by the epidemic
123136
uk_population <- epidemics::population(
124137
name = "UK",
125-
contact_matrix = contact_matrix,
138+
contact_matrix = contacts_byage_matrix,
126139
demography_vector = demography_vector,
127140
initial_conditions = initial_conditions
128141
)

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" "9711bc6201ff0d2fad7016a9e2402f2c" "site/built/contact-matrices.md" "2026-03-20"
8+
"episodes/simulating-transmission.Rmd" "6c36cc0947f2b782eac660b30b7572a5" "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" "ee2a68f0c593f53633226a38741a7a61" "site/built/modelling-interventions.md" "2026-03-20"
11+
"episodes/compare-interventions.Rmd" "43bfbe73569112710ca72808d18defc3" "site/built/compare-interventions.md" "2026-03-20"
12+
"episodes/vaccine-comparisons.Rmd" "b1b2e175887f763061a6c8bf900fe877" "site/built/vaccine-comparisons.md" "2026-03-20"
13+
"episodes/disease-burden.Rmd" "61be4094b1ecb9c2342f4bbd8e002e97" "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: 25 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,34 @@ 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)
94+
```
8995

96+
``` warning
97+
Warning: No reference provided.
98+
```
99+
100+
``` r
90101
# generate contact matrix
91-
cm_results <- socialmixr::contact_matrix(
92-
survey = survey_data,
102+
contacts_byage <- socialmixr::contact_matrix(
103+
survey = survey_load,
93104
countries = "United Kingdom",
94105
age_limits = c(0, 15, 65),
95-
symmetric = TRUE
106+
symmetric = TRUE,
107+
return_demography = TRUE
96108
)
97109

98110
# transpose contact matrix
99-
cm_matrix <- t(cm_results$matrix)
111+
contacts_byage_matrix <- t(contacts_byage$matrix)
100112

101113
# prepare the demography vector
102-
demography_vector <- cm_results$demography$population
103-
names(demography_vector) <- rownames(cm_matrix)
114+
demography_vector <- contacts_byage$demography$population
115+
names(demography_vector) <- rownames(contacts_byage_matrix)
104116

105117
# initial conditions: one in every 1 million is infected
106118
initial_i <- 1e-6
@@ -118,12 +130,12 @@ initial_conditions <- base::rbind(
118130
initial_conditions,
119131
initial_conditions
120132
)
121-
rownames(initial_conditions) <- rownames(cm_matrix)
133+
rownames(initial_conditions) <- rownames(contacts_byage_matrix)
122134

123135
# prepare the population to model as affected by the epidemic
124136
uk_population <- epidemics::population(
125137
name = "UK",
126-
contact_matrix = cm_matrix,
138+
contact_matrix = contacts_byage_matrix,
127139
demography_vector = demography_vector,
128140
initial_conditions = initial_conditions
129141
)
@@ -177,7 +189,7 @@ To include an intervention in our model we must create an `intervention` object.
177189

178190

179191
``` r
180-
rownames(cm_matrix)
192+
rownames(contacts_byage_matrix)
181193
```
182194

183195
``` output
@@ -419,8 +431,8 @@ Here we will assume all age groups are vaccinated at the same rate 0.01 and that
419431
# prepare a vaccination object
420432
vaccinate <- epidemics::vaccination(
421433
name = "vaccinate all",
422-
time_begin = matrix(40, nrow(cm_matrix)),
423-
time_end = matrix(40 + 150, nrow(cm_matrix)),
434+
time_begin = matrix(40, nrow(contacts_byage_matrix)),
435+
time_end = matrix(40 + 150, nrow(contacts_byage_matrix)),
424436
nu = matrix(c(0.01, 0.01, 0.01))
425437
)
426438
```

0 commit comments

Comments
 (0)