Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

R implementation of final size with susceptibility groups #42

Merged
merged 21 commits into from
Oct 3, 2022

Conversation

pratikunterwegs
Copy link
Member

@pratikunterwegs pratikunterwegs commented Sep 20, 2022

Adds functionality that allows final sizes in demographic groups with heterogeneous mixing and heterogeneous susceptibility. In short, demographic groups can be spread into discrete 'risk groups' (or susceptibility groups), each of which has a (potentially) different susceptibility to the epidemic. This fixes #41.

This PR adds the function final_size_grps which fixes #45, which has the following additional arguments:

  1. susceptibility: A matrix of susceptibility to the epidemic, where matrix element [i,j] represents the susceptibility of demographic group i in risk group j.
  2. p_susceptibility: A matrix representing the proportion of (or probability of) demographic group i being in risk group j.

Includes two internal functions, epi_spread and solve_final_size_iterative, which correspond to internal functions in #39. This code is essentially an R implementation of https://gitlab.com/epidemics-r/code_snippets/-/blob/feature/newton_solver/include/finalsize.hpp#L20 etc. These fix #44.

Also includes more tests for all three functions; fixes #46, and fixes #47. In general, input checking is restricted to top level functions exposed to the user, while unit tests mostly target internal functions (e.g. checking whether final size per group is in the range [0, 1]). Tests are based on https://gitlab.com/epidemics-r/code_snippets/-/blob/feature/newton_solver/tests/catch_final_size.cpp. Documentation added fixes #48.

@pratikunterwegs pratikunterwegs self-assigned this Sep 20, 2022
@pratikunterwegs pratikunterwegs added the New feature New feature or request label Sep 20, 2022
@pratikunterwegs pratikunterwegs marked this pull request as ready for review September 29, 2022 14:47
@pratikunterwegs
Copy link
Member Author

Hi @thimotei and @joshwlambert, I've requested that you review this PR. The idea is 1. to actually find things to fix and improve before merging into main, but also 2. to familiarise you with code review, and with the finalsize package.

See a little discussion on code review for PRs here: epiverse-trace/blueprints#6 (comment)

The PR references a number of issues, and the main comment above lays out a summary of the PR, and how it addresses the linked issues. Hopefully this is informative - if not, add a comment in the thread. If code in any of the files is unclear, add a comment on that chunk, which will open a conversation, to which I will reply.

If you can think of an easy fix or improvement, make a suggestion for an edit, which I can look at and potentially accept as is; this will make you co-authors on the PR.

If you have any questions about reviewing, you can ask me, @Bisaloo or @thibautjombart.

R/epi_spread.R Outdated Show resolved Hide resolved
R/iterative_solver.R Outdated Show resolved Hide resolved
R/iterative_solver.R Outdated Show resolved Hide resolved
@TimTaylor
Copy link

TimTaylor commented Sep 29, 2022

I've had a quick skim - will leave a more thorough review to others. Minor comments above. A couple of general points:

  • when porting the C++ code over to R, it's probably worth trying to make use of R's vectorisation as that will remove the need for some of the for loops and improve the performance.
  • I think some of the 1d matrices could actually be treated as dimensionless vectors. It looks like you convert them via as.vector in multiple places and I'm not sure if you ever use them as part of a matrix multiplication. If not it may simplify the code somewhat to just have them as vectors.

@pratikunterwegs
Copy link
Member Author

I've had a quick skim - will leave a more thorough review to others. Minor comments above. A couple of general points ...

Thanks @TimTaylor! Will work on the vectorisation and treating 1D matrices as vectors - as you noticed it's a direct translation from the Cpp - could def be something to be gained there.

#' data = 1, nrow = n_demo_grps, ncol = n_risk_grps
#' )
#' final_size_grps(
#' contact_matrix = r0 * c_matrix,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Comment on the readability: here c_matrix is named, could also potentially name d_vector, psusc and susc to make the mapping of age groups explicit between objects, more personal preference than anything.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So the c_matrix comes from socialmixr::polymod and carries over those naming conventions. I'm not sure why the d_vector, which also comes from the polymod data, is not named. Naming rows and columns for clarity is a good idea, but it's unclear where the naming should happen (in socialmixr? in finalsize?), and what should happen if none of the elements has names at all. I'm leaving this for now but it's worth perhaps opening an issue about it.

@joshwlambert
Copy link
Member

Only a couple of very minor points raised in the code. Overall, code is clean and works as expected. I have not fully analysed the details of the model but rather the code itself. Happy to approve. Nice work @pratikunterwegs.

@pratikunterwegs
Copy link
Member Author

Only a couple of very minor points raised in the code. Overall, code is clean and works as expected. I have not fully analysed the details of the model but rather the code itself. Happy to approve. Nice work @pratikunterwegs.

Thanks again @joshwlambert, and thanks for catching the superfluous test file and for noticing that pi overwrites the inbuilt value. I've changed that to something more informative as well. I'll keep this open until Monday morning, then merge it. Perhaps we can go over rebasing the feature/susc_grps branch on main after this merge, that should be a challenge. :)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
New feature New feature or request
3 participants