diff options
Diffstat (limited to 'R')
| -rw-r--r-- | R/ID.R | 13 | ||||
| -rw-r--r-- | R/IDEA.R | 13 | ||||
| -rw-r--r-- | R/WP.R | 29 | ||||
| -rw-r--r-- | R/WP_known.R | 7 | ||||
| -rw-r--r-- | R/WP_unknown.R | 14 | ||||
| -rw-r--r-- | R/computeLL.R | 9 | ||||
| -rw-r--r-- | R/seqB.R | 33 |
7 files changed, 58 insertions, 60 deletions
@@ -9,15 +9,14 @@ #' This method is based on an approximation of the SIR model, which is most valid at the beginning of an epidemic.
#' The method assumes that the mean of the serial distribution (sometimes called the serial interval) is known.
#' The final estimate can be quite sensitive to this value, so sensitivity testing is strongly recommended.
-#' Users should be careful about units of time (e.g. are counts observed daily or weekly?) when implementing.
+#' Users should be careful about units of time (e.g., are counts observed daily or weekly?) when implementing.
#'
-#' @param NT Vector of case counts
-#' @param mu Mean of the serial distribution (needs to match case counts in time units; for example, if case counts are
-#' weekly and the serial distribution has a mean of seven days, then \code{mu} should be set to one, if case
-#' counts are daily and the serial distribution has a mean of seven days, then \code{mu} should be set to seven)
+#' @param NT Vector of case counts.
+#' @param mu Mean of the serial distribution. This needs to match case counts in time units. For example, if case counts
+#' are weekly and the serial distribution has a mean of seven days, then \code{mu} should be set to one If case
+#' counts are daily and the serial distribution has a mean of seven days, then \code{mu} should be set to seven.
#'
-#' @return \code{ID} returns a list containing the following components: \code{Rhat} is the estimate of R0 and
-#' \code{inputs} is a list of the original input variables \code{NT, mu}.
+#' @return \code{ID} returns a single value, the estimate of R0.
#'
#' @examples
#' ## ===================================================== ##
@@ -10,15 +10,14 @@ #' This method is based on an approximation of the SIR model, which is most valid at the beginning of an epidemic.
#' The method assumes that the mean of the serial distribution (sometimes called the serial interval) is known.
#' The final estimate can be quite sensitive to this value, so sensitivity testing is strongly recommended.
-#' Users should be careful about units of time (e.g. are counts observed daily or weekly?) when implementing.
+#' Users should be careful about units of time (e.g., are counts observed daily or weekly?) when implementing.
#'
-#' @param NT Vector of case counts
-#' @param mu Mean of the serial distribution (needs to match case counts in time units; for example, if case counts are
-#' weekly and the serial distribution has a mean of seven days, then \code{mu} should be set to one, if case
-#' counts are daily and the serial distribution has a mean of seven days, then \code{mu} should be set to seven)
+#' @param NT Vector of case counts.
+#' @param mu Mean of the serial distribution. This needs to match case counts in time units. For example, if case counts
+#' are weekly and the serial distribution has a mean of seven days, then \code{mu} should be set to one. If case
+#' counts are daily and the serial distribution has a mean of seven days, then \code{mu} should be set to seven.
#'
-#' @return \code{IDEA} returns a list containing the following components: \code{Rhat} is the estimate of R0 and
-#' \code{inputs} is a list of the original input variables \code{NT, mu}.
+#' @return \code{IDEA} returns a single value, the estimate of R0.
#'
#' @examples
#' ## ===================================================== ##
@@ -19,7 +19,7 @@ source("WP_unknown.R") #' is strongly recommended. If the serial distribution is unknown (i.e., \code{mu} is \code{NA}), the #' likelihood function can be flat near the maximum, resulting in numerical instability of the optimizer. #' When \code{mu} is \code{NA}, the implementation takes considerably longer to run. Users should be careful -#' about units of time (e.g. are counts observed daily or weekly?) when implementing. +#' about units of time (e.g., are counts observed daily or weekly?) when implementing. #' #' The model developed in White and Pagano (2008) is discrete, and hence the serial distribution is finite #' discrete. In our implementation, the input value \code{mu} is that of a continuous distribution. The @@ -28,27 +28,26 @@ source("WP_unknown.R") #' user notices that the input \code{mu} and output mean of \code{SD} are different, this is to be expected, #' and is caused by the discretization. #' -#' @param NT Vector of case counts +#' @param NT Vector of case counts. #' @param mu Mean of the serial distribution (needs to match case counts in time units; for example, if case #' counts are weekly and the serial distribution has a mean of seven days, then \code{mu} should be #' set to one). The default value of \code{mu} is set to \code{NA}. -#' @param search List of default values for the grid search algorithm; the list includes three elements: the -#' first is \code{B} which is the length of the grid in one dimension, the second is -#' \code{scale.max} which is the largest possible value of the scale parameter, and the third is -#' \code{shape.max} which is the largest possible value of the shape parameter; defaults to +#' @param search List of default values for the grid search algorithm. The list includes three elements: the +#' first is \code{B}, which is the length of the grid in one dimension; the second is +#' \code{scale.max}, which is the largest possible value of the scale parameter; and the third +#' is \code{shape.max}, which is the largest possible value of the shape parameter. Defaults to #' \code{B=100, scale.max=10, shape.max=10}. For both shape and scale, the smallest possible -#' value is 1/\code{B}. +#' value is 1/\code{B}. #' @param tol Cutoff value for cumulative distribution function of the pre-discretization gamma serial -#' distribution, defaults to 0.999 (i.e. in the discretization, the maximum is chosen such that the +#' distribution. Defaults to 0.999 (i.e. in the discretization, the maximum is chosen such that the #' original gamma distribution has cumulative probability of no more than 0.999 at this maximum). #' -#' @return WP returns a list containing the following components: \code{Rhat} is the estimate of R0, \code{SD} -#' is either the discretized serial distribution (if \code{mu} is not \code{NA}) or the estimated -#' discretized serial distribution (if \code{mu} is \code{NA}), and \code{inputs} is a list of the -#' original input variables \code{NT, mu, method, search, tol}. The list also returns the variable -#' \code{check}, which is equal to the number of non-unique maximum likelihood estimators. The serial -#' distribution \code{SD} is returned as a list made up of \code{supp} the support of the distribution -#' and \code{pmf} the probability mass function. +#' @return \code{WP} returns a list containing the following components: \code{Rhat} is the estimate of R0, +#' and \code{SD} is either the discretized serial distribution (if \code{mu} is not \code{NA}), or the +#' estimated discretized serial distribution (if \code{mu} is \code{NA}). The list also returns the +#' variable \code{check}, which is equal to the number of non-unique maximum likelihood estimators. +#' The serial distribution \code{SD} is returned as a list made up of \code{supp} (the support of +#' the distribution) and \code{pmf} (the probability mass function). #' #' @examples #' ## ===================================================== ## diff --git a/R/WP_known.R b/R/WP_known.R index 563b7ae..2c54bce 100644 --- a/R/WP_known.R +++ b/R/WP_known.R @@ -3,9 +3,10 @@ #' This is a background/internal function called by \code{WP}. It computes the maximum #' likelihood estimator of R0 assuming that the serial distribution is known and finite discrete. #' -#' @param NT vector of case counts -#' @param p discretized version of the serial distribution -#' @return The function returns \code{Rhat}, the maximum likelihood estimator of R0. +#' @param NT Vector of case counts. +#' @param p Discretized version of the serial distribution. +#' +#' @return The function returns the maximum likelihood estimator of R0. #' #' @export WP_known <- function(NT, p) { diff --git a/R/WP_unknown.R b/R/WP_unknown.R index 569a880..542f0f7 100644 --- a/R/WP_unknown.R +++ b/R/WP_unknown.R @@ -8,11 +8,11 @@ source("WP_known.R") #' The function then implements a simple grid search algorithm to obtain the maximum likelihood estimator #' of R0 as well as the gamma parameters. #' -#' @param NT vector of case counts -#' @param B length of grid for shape and scale (grid search parameter) -#' @param shape.max maximum shape value (grid search parameter) -#' @param scale.max maximum scale value (grid search parameter) -#' @param tol cutoff value for cumulative distribution function of the serial distribution, defaults to 0.999 +#' @param NT Vector of case counts. +#' @param B Length of grid for shape and scale (grid search parameter). +#' @param shape.max Maximum shape value (grid \code{search} parameter). +#' @param scale.max Maximum scale value (grid \code{search} parameter). +#' @param tol cutoff value for cumulative distribution function of the serial distribution (defaults to 0.999). #' #' @return The function returns \code{Rhat}, the maximum likelihood estimator of R0, as well as the maximum #' likelihood estimator of the discretized serial distribution given by \code{p} (the probability mass @@ -20,8 +20,8 @@ source("WP_known.R") #' The function also returns \code{resLL} (all values of the log-likelihood) at \code{shape} (grid for #' shape parameter) and at \code{scale} (grid for scale parameter), as well as \code{resR0} (the full #' vector of maximum likelihood estimators), \code{JJ} (the locations for the likelihood for these), and -#' \code{J0} (the location for the maximum likelihood estimator \code{Rhat}). If \code{JJ} and \code{J0} -#' are not the same, this means that the maximum likelihood estimator is not unique. +#' \code{J0} (the location for the maximum likelihood estimator \code{Rhat}). If \code{JJ} and \code{J0} +#' are not the same, this means that the maximum likelihood estimator is not unique. #' #' @export WP_unknown <- function(NT, B=100, shape.max=10, scale.max=10, tol=0.999) { diff --git a/R/computeLL.R b/R/computeLL.R index 6acfc78..9047777 100644 --- a/R/computeLL.R +++ b/R/computeLL.R @@ -2,10 +2,11 @@ #' #' This is a background/internal function called by \code{WP}. It computes the log-likelihood. #' -#' @param NT vector of case counts -#' @param p discretized version of the serial distribution -#' @param R0 basic reproductive ratio -#' @return The function returns the variable \code{LL} which is the log-likelihood at the input variables and parameters. +#' @param NT Vector of case counts. +#' @param p Discretized version of the serial distribution. +#' @param R0 Basic reproductive ratio. +#' +#' @return This function returns the log-likelihood at the input variables and parameters. #' #' @export computeLL <- function(p, NT, R0) { @@ -7,43 +7,42 @@ #' The distribution of R0 is then updated sequentially, with one update for each new case count observation. #' The final estimate of R0 is \code{Rhat}, the mean of the (last) posterior distribution. #' The prior distribution is the initial belief of the distribution of R0; which in this implementation is the uninformative uniform -#' distribution with values between zero and \code{kappa}. Users can change the value of kappa only (ie. the prior distribution +#' distribution with values between zero and \code{kappa}. Users can change the value of /code{kappa} only (i.e., the prior distribution #' cannot be changed from the uniform). As more case counts are observed, the influence of the prior distribution should lessen on #' the final estimate \code{Rhat}. #' #' This method is based on an approximation of the SIR model, which is most valid at the beginning of an epidemic. The method assumes #' that the mean of the serial distribution (sometimes called the serial interval) is known. The final estimate can be quite sensitive -#' to this value, so sensitivity testing is strongly recommended. Users should be careful about units of time (e.g. are counts observed +#' to this value, so sensitivity testing is strongly recommended. Users should be careful about units of time (e.g., are counts observed #' daily or weekly?) when implementing. #' #' Our code has been modified to provide an estimate even if case counts equal to zero are present in some time intervals. This is done #' by grouping the counts over such periods of time. Without grouping, and in the presence of zero counts, no estimate can be provided. #' -#' @param NT Vector of case counts -#' @param mu Mean of the serial distribution (needs to match case counts in time units; for example, if case counts are -#' weekly and the serial distribution has a mean of seven days, then \code{mu} should be set to one, if case -#' counts are daily and the serial distribution has a mean of seven days, then \code{mu} should be set to seven) -#' @param kappa Largest possible value of uniform prior, defaults to 20. This describes the prior belief on ranges of R0, -#' so should be set to a higher value if R0 is believed to be larger. +#' @param NT Vector of case counts. +#' @param mu Mean of the serial distribution. This needs to match case counts in time units. For example, if case counts +#' are weekly and the serial distribution has a mean of seven days, then \code{mu} should be set to one. If case +#' counts are daily and the serial distribution has a mean of seven days, then \code{mu} should be set to seven. +#' @param kappa Largest possible value of uniform prior (defaults to 20). This describes the prior belief on ranges of R0, +#' and should be set to a higher value if R0 is believed to be larger. #' -#' @return secB returns a list containing the following components: \code{Rhat} is the estimate of R0 (the posterior mean), -#' \code{posterior} is the posterior distribution of R0 from which alternate estimates can be obtained (see examples), -#' \code{group} is an indicator variable (if \code{group=TRUE}, zero values of NT were input and grouping was done to -#' obtain \code{Rhat}), and \code{inputs} is a list of the original input variables \code{NT, gamma, kappa}. The variable -#' \code{posterior} is returned as a list made up of \code{supp} the support of the distribution and \code{pmf} the -#' probability mass function. +#' @return \code{secB} returns a list containing the following components: \code{Rhat} is the estimate of R0 (the posterior mean), +#' \code{posterior} is the posterior distribution of R0 from which alternate estimates can be obtained (see examples), +#' and \code{group} is an indicator variable (if \code{group=TRUE}, zero values of NT were input and grouping was done +#' to obtain \code{Rhat}). The variable \code{posterior} is returned as a list made up of \code{supp} (the support of +#' the distribution) and \code{pmf} (the probability mass function). #' #' @examples #' ## ===================================================== ## #' ## Illustrate on weekly data ## #' ## ===================================================== ## #' -#' NT <- c(1, 4, 10, 5, 3, 4, 19, 3, 3, 14, 4) +#' NT <- c(1, 4, 10, 5, 3, 4, 19, 3, 3, 14, 4) #' ## obtain Rhat when serial distribution has mean of five days -#' res1 <- seqB(NT=NT, mu=5/7) +#' res1 <- seqB(NT=NT, mu=5/7) #' res1$Rhat #' ## obtain Rhat when serial distribution has mean of three days -#' res2 <- seqB(NT=NT, mu=3/7) +#' res2 <- seqB(NT=NT, mu=3/7) #' res2$Rhat #' #' ## ============================================================= ## |