/
parameters_vaccine.R
169 lines (139 loc) · 6.68 KB
/
parameters_vaccine.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
# --------------------------------------------------------------------------------
# infection process for vaccination model (multiple doses, no types)
# Sean L. Wu (slwood89@gmail.com)
# July 2021
# --------------------------------------------------------------------------------
#' @title Get vaccine efficacy and Ab titre parameters
#' @param vaccine which vaccine? should be one of: "Pfizer", "AstraZeneca", "Sinovac", "Moderna" (controls the mean titre associated with each dose)
#' @param max_dose maximum number of doses
#' @param correlated are doses correlated?
#' @param hl_s Half life of antibody decay - short
#' @param hl_l Half life of antibody decay - long
#' @param period_s length of the initial decay rate (shortest half-life)
#' @param t_period_l Time point at which to have switched to longest half-life
#' @param ab_50 titre relative to convalescent required to provide 50% protection from infection, on linear scale
#' @param ab_50_severe titre relative to convalescent required to provide 50% protection from severe disease, on linear scale
#' @param std10 Pooled standard deviation of antibody level on log10 data
#' @param k shape parameter of efficacy curve
#' @param nt_transmission_factor used by [safir::vaccine_efficacy_transmission] to compute the effect of antibody titre on onward transmission
#' @param nt_efficacy_transmission a logical variable, specify if antibody titre should affect onward transmission
#' @param max_ab maximum allowable antibody titre draw (on natural log scale)
#' @param mu_ab_list a data.frame
#' @description Get parameters for vaccine efficacy and antibody titre decay rate.
#' @export
get_vaccine_ab_titre_parameters <- function(
vaccine, max_dose = 2, correlated = FALSE,
hl_s = 108, hl_l = 3650,
period_s = 250, t_period_l = 365,
ab_50 = 0.2,
ab_50_severe = 0.03,
std10 = 0.44,
k = 2.94,
nt_transmission_factor = 12,
nt_efficacy_transmission = FALSE,
max_ab = 8,
mu_ab_list = data.frame(name = c("Pfizer", "AstraZeneca", "Sinovac", "Moderna"),
mu_ab_d1 = c(13/94, 1/59, 28/164, ((185+273)/2)/321),
mu_ab_d2 = c(223/94, 32/59, 28/164, 654/158))
) {
stopifnot(is.logical(nt_efficacy_transmission))
stopifnot(max_dose %in% 1:(ncol(mu_ab_list) - 1))
stopifnot(is.logical(correlated))
stopifnot(vaccine %in% mu_ab_list[, "name"])
stopifnot(is.finite(nt_transmission_factor))
stopifnot(nt_transmission_factor > 0)
stopifnot(is.finite(max_ab))
stopifnot(max_ab > 0)
stopifnot(all(is.finite(c(hl_s, hl_l, period_s, t_period_l, ab_50, ab_50_severe, std10, k))))
stopifnot(all(c(hl_s, hl_l, period_s, t_period_l, ab_50, ab_50_severe, std10, k) > 0))
time_to_decay <- t_period_l - period_s # time in days to reach longest half-life
dr_s <- -log(2)/hl_s # Corresponding decay rate in days for half life above
dr_l <- -log(2)/hl_l
dr_vec <- c(
rep(dr_s, period_s),
seq(dr_s, dr_l, length.out = time_to_decay),
dr_l
)
mu_ab <- as.numeric(mu_ab_list[mu_ab_list$name == vaccine, -1])
stopifnot(length(mu_ab) >= max_dose)
stopifnot(all(is.finite(mu_ab)))
parameters <- list(
# necessary
dr_vec = dr_vec,
mu_ab = mu_ab,
std10 = std10,
ab_50 = ab_50,
ab_50_severe = ab_50_severe,
k = k,
nt_transmission_factor = nt_transmission_factor,
nt_efficacy_transmission = nt_efficacy_transmission,
max_ab = max_ab,
correlated = correlated,
# other
hl_s = hl_s,
hl_l = hl_l,
period_s = period_s,
t_period_l = t_period_l
)
return(parameters)
}
#' @title Combine and verify vaccine parameters
#' @note If modeling a single dose, `dose_period` must be a vector of length 1 and
#' `next_dose_priority_matrix` may be set to `NULL`.
#' @param safir_parameters a list from \code{\link{get_parameters}}
#' @param vaccine_ab_parameters a list from \code{\link{get_vaccine_ab_titre_parameters}}
#' @param vaccine_set a vector giving the number of doses available each day (not each timestep)
#' @param dose_period a vector giving the minimum delay between doses
#' @param strategy_matrix a vaccine strategy matrix from \code{\link[nimue]{strategy_matrix}}
#' or a list of strategy matrices, one for each vaccine dose
#' @param next_dose_priority_matrix a binary matrix giving age groups prioritized for next dose;
#' it should have one fewer row than the number of doses being given, because on the
#' final allocation phase there will be no future dose to prioritize
#' @description Combine parameters for simulation and verify for correctness.
#' @export
make_vaccine_parameters <- function(safir_parameters, vaccine_ab_parameters, vaccine_set, dose_period, strategy_matrix, next_dose_priority_matrix) {
parameters <- safir_parameters
vaccine_doses <- length(dose_period)
stopifnot(vaccine_doses >= 1)
if (vaccine_doses > 1) {
stopifnot(all(is.finite(dose_period[2:vaccine_doses])))
stopifnot(all(dose_period[2:vaccine_doses] >= 0))
}
parameters$N_phase <- vaccine_doses
parameters$dose_period <- dose_period
stopifnot(length(vaccine_set) == parameters$time_period)
stopifnot(all(is.finite(vaccine_set)))
stopifnot(all(vaccine_set >= 0))
parameters$vaccine_set <- floor(vaccine_set)
storage.mode(next_dose_priority_matrix) <- "integer"
stopifnot(nrow(next_dose_priority_matrix) == vaccine_doses - 1)
# one strategy matrix for all dose phases (input matrix)
if (inherits(strategy_matrix, 'matrix')) {
N_age <- ncol(strategy_matrix)
strategy_matrix_array <- replicate(n = vaccine_doses, expr = strategy_matrix, simplify = FALSE)
# one strategy matrix for each dose phase (input list)
} else if(inherits(strategy_matrix, 'list')) {
stopifnot(length(strategy_matrix) == vaccine_doses)
N_age <- ncol(strategy_matrix[[1]])
strategy_matrix_array <- strategy_matrix
# bad input
} else {
stop("invalid object passed for strategy_matrix")
}
stopifnot(ncol(next_dose_priority_matrix) == N_age)
stopifnot(parameters$N_age == N_age)
parameters$vaccine_coverage_mat <- strategy_matrix_array
parameters$next_dose_priority <- next_dose_priority_matrix
# attach necessary parameters from vaccine pars
parameters$dr_vec <- vaccine_ab_parameters$dr_vec
parameters$mu_ab <- vaccine_ab_parameters$mu_ab
parameters$std10 <- vaccine_ab_parameters$std10
parameters$ab_50 <- vaccine_ab_parameters$ab_50
parameters$ab_50_severe <- vaccine_ab_parameters$ab_50_severe
parameters$k <- vaccine_ab_parameters$k
parameters$nt_transmission_factor <- vaccine_ab_parameters$nt_transmission_factor
parameters$nt_efficacy_transmission <- vaccine_ab_parameters$nt_efficacy_transmission
parameters$max_ab <- vaccine_ab_parameters$max_ab
parameters$correlated <- vaccine_ab_parameters$correlated
return(parameters)
}