---
title: "Modelling multiple contact settings"
bibliography: resources/references.json
link-citations: true
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Modelling multiple contact settings}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

This vignette shows how to pass community contact matrices disaggregated into multiple settings to the Daedalus model.

Daedalus allows users to apply closures or any arbitrary scaling to each setting separately, rather than uniformly scaling contacts across settings.
An example of setting-specific scaling is scaling school contacts for school-age children, while leaving contacts in other settings intact.

```{r setup}
# load package
library(daedalus)
library(socialmixr) # for POLYMOD data
```

## Default contact settings

Daedalus countries have two contact settings by default as of version 0.3.6:

1. Community or 'total': This is the set of age-specific contacts that take place in the community, outside the workplace. In daedalus < v0.3.5, this would have been the element referred to as `x$contact_matrix`;

2. Workplace: This is a matrix giving the contacts within workplaces (so between workers in an economic sector; forming the diagonal of the matrix), and from non-working consumers to workers in each economic sector (considers only contacts from consumers to workers and not vice versa).

These are displayed in the country object when printed to screen.

```{r default_country}
daedalus_country("GBR")
```

## Basic example using POLYMOD

This example using POLYMOD data [@mossong2008] shows how to assign contacts from two different settings to a `<daedalus_country>` object.

Once this assignment is done, helper functions such as `prepare_parameters()` take care of expanding and scaling the contact matrices when the country object is passed to `daedalus()`.

**Note that** running `daedalus("XYZ", ...)` cannot be used to model multiple settings, and the `country` argument must be a `<daedalus_country>` if multiple settings are needed.

We can get POLYMOD data [@mossong2008] for the U.K. from the [R package _socialmixr_](https://cran.r-project.org/package=socialmixr).
Here in this basic example, we consider only two settings, 'home' and 'school'.

```{r get_polymod}
# load polymod data and get home and school contacts
data(polymod)

# get matrices with Daedalus age bins
cm_uk_home <- contact_matrix(
  polymod,
  "GB",
  age_limits = c(0, 5, 18, 65),
  filter = list(cnt_home = 1)
)$matrix

cm_uk_school <- contact_matrix(
  polymod,
  "GB",
  age_limits = c(0, 5, 18, 65),
  filter = list(cnt_school = 1)
)$matrix
```

Create a new country object, passing the new settings to the `contact_matrix` argument as a named list.
Having names helps when examining the object later.

```{r multi_setting}
x <- daedalus_country(
  "GB",
  contact_matrix = list(
    home = cm_uk_home,
    school = cm_uk_school
  )
)

# print x to see output
x
```

Pass the created country to [daedalus()] as usual to run the epi-econ model taking into account contacts in each setting.

```{r run_daedalus}
daedalus(x, "sars_cov_1", time_end = 100)
```

## Passing setting-specific contact scaling

Users can pass setting-specific contact scaling by using the `<daedalus_npi>` class.

Scaling must take the form of an $M \times N$ matrix, where $M$ is the number of age and economic strata (49 for Daedalus), and $N$ is the number of settings.
For a default `<daedalus_country>`, this would be a $49 \times 2$ matrix, due to the two contact settings, 'community' and 'workplace'.

The example below shows a matrix that does not scale community contacts, and which scales all workplace contacts to 30% of the total.
The underlying model will scale the within-sector contacts by the square of the scaling coefficient (representing scaling of both incoming and outgoing contacts), while the consumer-to-worker contacts will be scaled by the coefficient alone, as we are not simulating a reduction in the number of consumers visiting workplaces.

```{r scaling_settings}
openness <- cbind(
  rep(1, 49), # no scaling on community
  c(rep(1, 4), rep(0.3, 45)) # scale 45 economic strata
)
```

The following construction is simpler because it assumes the same scaling for age-specific contacts in workplaces as for worker contacts, and would not cause issues as age-specific workplace contacts are 0.
Nonetheless, we advise using the construction above as it makes the intention more explicit.

```{r scaling_settings_2}
openness <- cbind(
  rep(1, 49), # no scaling on community
  rep(0.3, 49) # scale all 49 age-and-economic strata
)
```

A third alternative is to take advantage of the default behaviour of the `<daedalus_npi>` class constructor, which is to assume that the scaling of community contacts is the simple mean of the scaling of workplace contacts, and pass a single 45-element workplace contacts scaling vector.
This also applies to [daedalus_timed_npi()].

```{r scaling_settings_3}
openness <- rep(0.3, 45) # for 45 economic sectors
```

Create a `<daedalus_npi>`, passing the constructed matrix `openness`, and use it with `daedalus()` as needed.

```{r custom_npi}
# a 90-day timed NPI operating on GB; note timed_npi needs a list of openness
daedalus_timed_npi(30, 90, list(openness), "GB")
```

Openness values for pre-defined NPIs are also stored in this way.

**Note that** pre-defined NPIs include community contacts scaling by default, as shown in the code snippet below.
The scaling of community contacts is approximately 0.71, which is the mean of the scaling of workplace contacts shown below.

```{r predef_npis}
x <- daedalus_npi("elimination", "GB", "sars_cov_1")
x

x$parameters$openness$openness
```

## Future developments

This functionality is still in development, and features being considered include applying different behavioural modifiers to each setting to simulate endogenous risk-perception related modification in contacts across settings.

## References