Update DESCRIPTION

[Rnaught] / R / IDEA.R

diff --git a/R/IDEA.R b/R/IDEA.R
index 9ca1da67e9721d1116ce40828231827ea3775e5e..854acd7af108903da248048d34ee6ae20de3a81d 100644 (file)
--- a/R/IDEA.R
+++ b/R/IDEA.R
@@ -1,63 +1,58 @@
 #' IDEA method
 #'
 #' IDEA method
 #'

-#' This function implements a least squares estimation method of R0 due to Fisman et al. (PloS One, 2013).    See details for implementation notes.
+#' This function implements a least squares estimation method of R0 due to Fisman et al. (PloS One, 2013).
+#' See details for implementation notes.

 #'
 #'

-#'This method is closely related to that implemented in \code{ID}. The method is based on an incidence decay model.   The estimate of R0 is the value which minimizes the sum of squares between observed case counts and cases counts 'expected' under the model.
+#' This method is closely related to that implemented in \code{ID}. The method is based on an incidence decay model.
+#' The estimate of R0 is the value which minimizes the sum of squares between observed case counts and cases counts
+#' expected under the model.

 #'
 #'

-#' 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 recommended. Users should be careful about units of time (e.g. are counts observed daily or weekly?) when implementing.  
+#' 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.

 #'
 #'

-#' @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)
+#' @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{IDEA} returns a single value, the estimate of R0.

 #'
 #' @examples
 #'
 #' @examples

-#' 

 #' ## ===================================================== ##
 #' ## Illustrate on weekly data                             ##
 #' ## ===================================================== ##
 #'
 #' ## ===================================================== ##
 #' ## 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
 #' ## obtain Rhat when serial distribution has mean of five days

-#' res1 <- IDEA(NT=NT, mu=5/7) 
-#' res1$Rhat
+#' IDEA(NT=NT, mu=5/7)

 #' ## obtain Rhat when serial distribution has mean of three days
 #' ## obtain Rhat when serial distribution has mean of three days

-#' res2        <- IDEA(NT=NT, mu=3/7)  
-#' res2$Rhat
+#' IDEA(NT=NT, mu=3/7)

 #'
 #' ## ========================================================= ##
 #' ## Compute Rhat using only the first five weeks of data      ##
 #' ## ========================================================= ##
 #'
 #'
 #' ## ========================================================= ##
 #' ## Compute Rhat using only the first five weeks of data      ##
 #' ## ========================================================= ##
 #'

-#' 
-#' res3 <- IDEA(NT=NT[1:5], mu=5/7)            # serial distribution has mean of five days
-#' res3$Rhat
-#' @export
+#' IDEA(NT=NT[1:5], mu=5/7) # serial distribution has mean of five days

 #'
 #'

-\r
-\r
-IDEA <- function(NT, mu){\r
+#' @export
+IDEA <- function(NT, mu) {
+    if (length(NT) < 2)
+        print("Warning: length of NT should be at least two.")
+    else {
+        NT <- as.numeric(NT)
+        TT <- length(NT)
+        s <- (1:TT) / mu
+
+        y1 <- log(NT) / s
+        y2 <- s^2
+        y3 <- log(NT)

 
 

-       if(length(NT)<2) {
-               print("Warning: length of NT should be at least two.")
-       }
-       else{
-         NT            <-      as.numeric(NT)
-         TT            <-      length(NT)
-         s             <-      (1:TT)/mu
-         
-         y1            <-      log(NT)/s
-         y2            <-      s^2
-         y3            <-      log(NT)
-#        IDEA1         <-  cumsum(y2)*cumsum(y1)-cumsum(s)*cumsum(y3)
-#        IDEA2         <-      (1:TT)*cumsum(y2)-(cumsum(s))^2 
-#        IDEA          <-      exp(IDEA1/IDEA2)
-#        Rhat          <-      tail(IDEA,1)
-         IDEA1         <-  sum(y2)*sum(y1)-sum(s)*sum(y3)
-         IDEA2         <-      TT*sum(y2)-(sum(s))^2   
-         IDEA          <-      exp(IDEA1/IDEA2)
-         
+        IDEA1 <- sum(y2) * sum(y1) - sum(s) * sum(y3)
+        IDEA2 <- TT * sum(y2) - sum(s)^2
+        IDEA <- exp(IDEA1 / IDEA2)

 
 

-         return(list(Rhat=IDEA, inputs=list(NT=NT, mu=mu)))
-         }\r
-}\r
+        return(IDEA)
+    }
+}