batch update

STAPM · Sep 17, 2021 · 4d8376d · 4d8376d
1 parent 38ec9e8
commit 4d8376d
Show file tree

Hide file tree

Showing 21 changed files with 2,465 additions and 668 deletions.
diff --git a/NAMESPACE b/NAMESPACE
@@ -13,15 +13,16 @@ export(RRtob)
 export(TobAlcInt)
 export(TobLags)
 export(WArisk_acute)
+export(intervalprob)
 export(subgroupRisk)
 import(data.table)
 importFrom(Rdpack,reprompt)
 importFrom(data.table,":=")
 importFrom(data.table,.N)
 importFrom(data.table,.SD)
-importFrom(data.table,CJ)
 importFrom(data.table,copy)
 importFrom(data.table,data.table)
+importFrom(data.table,fifelse)
 importFrom(data.table,rbindlist)
 importFrom(data.table,setDT)
 importFrom(data.table,setnames)

diff --git a/R/PAFcalc.R b/R/PAFcalc.R
@@ -1,12 +1,15 @@
 
-#' Calculate Population Attributable Fractions \lifecycle{maturing}
+#' Calculate Population Attributable Fractions
 #' 
 #' Uses \code{RRFunc()} and \code{subgroupRisk()} to 
 #' calculate population attributable fractions 
 #' based on the survey data provided.
 #' 
 #' 
 #' @param data Data table of individual characteristics
+#' @param rrdata Optional - data table containing individual tobacco and alcohol consumption characteristics 
+#' with relative risks of disease already assigned. This could be useful for increasing efficiency - saving computer processing time. 
+#' Defaults to NULL.
 #' @param substance Whether to compute relative risks for just alcohol ("alc"),
 #' just tobacco ("tob") or joint risks for tobacco and alcohol ("tobalc").
 #' @param tob_include_risk_in_former_smokers Logical - whether the residual risks of smoking in former smokers
@@ -49,6 +52,7 @@
 #' 
 PAFcalc <- function(
   data,
+  rrdata = NULL,
   substance,
   tob_include_risk_in_former_smokers = TRUE,
   alc_protective = TRUE,
@@ -66,54 +70,65 @@ PAFcalc <- function(
 
   years <- min(data$year):max(data$year)
 
-  cat("Assigning relative risks\n")
-
-  for(y in years) {
-
-    #y <- years[1]
-
-    cat("\t", y, "\n")
+  if(is.null(rrdata)) {
 
-    # Add the relative risks to the data
-    data_rr <- tobalcepi::RRFunc(
-      data = data[year == y],
-      substance = substance,
-      tob_diseases = tobalcepi::tob_disease_names,
-      tob_include_risk_in_former_smokers = tob_include_risk_in_former_smokers,
-      alc_diseases = tobalcepi::alc_disease_names,
-      alc_mort_and_morb = c(
-        "Ischaemic_heart_disease", 
-        "LiverCirrhosis", 
-        "Haemorrhagic_Stroke",
-        "Ischaemic_Stroke"),
-      alc_risk_lags = FALSE,
-      alc_protective = alc_protective,
-      alc_wholly_chronic_thresholds = alc_wholly_chronic_thresholds,
-      alc_wholly_acute_thresholds = alc_wholly_acute_thresholds,
-      grams_ethanol_per_unit = grams_ethanol_per_unit,
-      show_progress = TRUE,
-      within_model = within_model,
-      tobalc_include_int = tobalc_include_int)
+    cat("Assigning relative risks\n")
 
-    if(y == years[1]) {
+    for(y in years) {
 
-      data_rr_comb <- copy(data_rr)
+      #y <- years[1]
 
-    } else {
+      cat("\t", y, "\n")
 
-      data_rr_comb <- rbindlist(list(data_rr_comb, copy(data_rr)), use.names = T)
+      # Add the relative risks to the data
+      data_rr <- tobalcepi::RRFunc(
+        data = data[year == y],
+        substance = substance,
+        tob_diseases = tobalcepi::tob_disease_names,
+        tob_include_risk_in_former_smokers = tob_include_risk_in_former_smokers,
+        alc_diseases = tobalcepi::alc_disease_names,
+        alc_mort_and_morb = c(
+          "Ischaemic_heart_disease", 
+          "LiverCirrhosis", 
+          "Haemorrhagic_Stroke",
+          "Ischaemic_Stroke"),
+        alc_risk_lags = FALSE,
+        alc_protective = alc_protective,
+        alc_wholly_chronic_thresholds = alc_wholly_chronic_thresholds,
+        alc_wholly_acute_thresholds = alc_wholly_acute_thresholds,
+        grams_ethanol_per_unit = grams_ethanol_per_unit,
+        show_progress = TRUE,
+        within_model = within_model,
+        tobalc_include_int = tobalc_include_int)
+
+      if(y == years[1]) {
+
+        data_rr_comb <- copy(data_rr)
+
+      } else {
+
+        data_rr_comb <- rbindlist(list(data_rr_comb, copy(data_rr)), use.names = T)
+
+      }
 
     }
 
-  }
 
+
+
+  } else {
+
+    data_rr_comb <- copy(rrdata)
+
+  }
 
   # If need morbidity relative risks
 
   if(mort_or_morb == "morb") {
     data_rr_comb <- data_rr_comb[ , c("Ischaemic_heart_disease", "LiverCirrhosis", "Haemorrhagic_Stroke", "Ischaemic_Stroke") := NULL]
     setnames(data_rr_comb, paste0(c("Ischaemic_heart_disease", "LiverCirrhosis", "Haemorrhagic_Stroke", "Ischaemic_Stroke"), "_morb"), c("Ischaemic_heart_disease", "LiverCirrhosis", "Haemorrhagic_Stroke", "Ischaemic_Stroke"))
   }
+
 
 
   # Calculate PAFs

diff --git a/R/PArisk.R b/R/PArisk.R
@@ -1,6 +1,6 @@
 
 
-#' Relative risks for alcohol-related injuries \lifecycle{stable}
+#' Relative risks for alcohol-related injuries
 #'
 #' Uses the 'new' binge model methods to calculate a relative risk
 #' for each individual for experiencing each cause during one year.
@@ -13,10 +13,8 @@
 #'  and the discussion paper by Hill-McManus 2014. The relative risks for alcohol-related injuries
 #'  are taken from Cherpitel et al 2015.
 #'
-#' @param SODMean Numeric vector - the average amount that each individual is expected to
-#' drink on a single drinking occasion.
-#' @param SODSDV Numeric vector - the standard deviation of the amount that each individual is expected to
-#' drink on a single drinking occasion.
+#' @param interval_prob_vec Column of vectors - the probabilities that each individual 
+#' drinks each amount of grams of ethanol (1:600) on a single drinking occasion.
 #' @param SODFreq Numeric vector - the expected number of drinking occasions that
 #' each individual has each week.
 #' @param Weight Numeric vector - each individual's body weight in kg.
@@ -33,7 +31,7 @@
 #'
 #' @return Returns a numeric vector of each individual's relative risk of the acute consequences of drinking.
 #' 
-#' @importFrom data.table := setDT setnames CJ
+#' @importFrom data.table := setDT setnames
 #' 
 #' @export
 #'
@@ -137,8 +135,7 @@
 #' }
 #'
 PArisk <- function(
-  SODMean = NULL,
-  SODSDV = NULL,
+  interval_prob_vec = NULL,
   SODFreq = NULL,
   Weight = NULL,
   Widmark_r = NULL,
@@ -149,7 +146,14 @@ PArisk <- function(
   getcurve = FALSE
 ) {
 
-  kn <- 600
+
+  # The max number of grams of ethanol per day drunk on a single drinking ocassion
+  # 600 is v large
+  # but this can influence the probability density
+  # but it makes things slow to have a large number
+  # try scaling back by 10%
+  kn <- 600 * 0.9
+
 
   grams_ethanol <- 1:kn
 
@@ -181,55 +185,8 @@ PArisk <- function(
     # Convert to hours
     Duration_h <- Duration_m / 60
 
-    #######################
-    # Calculate the cumulative probability distribution of each amount of alcohol (1 to 100 g) being drunk on an occasion
-    # x <- stats::pnorm(
-    #   grams_ethanol,
-    #   SODMean * grams_ethanol_per_unit, # mean
-    #   SODSDV * grams_ethanol_per_unit # variance
-    # )
-
-    # grams_ethanol <- 1:600
-    # SODMean <- 4
-    # SODSDV <- 2
-    # grams_ethanol_per_unit <- 8
-    # lb <- bench::mark(
-    # x <- t(sapply(grams_ethanol,
-    #               stats::pnorm, 
-    #               mean = SODMean * grams_ethanol_per_unit, # mean
-    #               sd = SODSDV * grams_ethanol_per_unit # variance
-    # ))
-    # ,
-    # 
-
-    x <- t(vapply(X = grams_ethanol,
-                  FUN = stats::pnorm, 
-                  FUN.VALUE = numeric(length(SODMean)),
-                  mean = SODMean * grams_ethanol_per_unit, # mean
-                  sd = SODSDV * grams_ethanol_per_unit # variance
-    ))
-
-    # 
-    # )
-    # 
-    # lb
-
-    #######################
-    # Convert from the cumulative distribution to the
-    # probability that each level of alcohol is consumed on a drinking occasion
-    #interval_prob <- x - c(0, x[1:(length(x) - 1)])
-    #interval_prob <- diff(x)
-    interval_prob <- apply(x, 2, diff)
-
-    #interval_prob <- interval_prob / sum(interval_prob)
-
-    # NOT SURE IF THE LINE BELOW IS NEEDED
-    # the code makes values sum to 1
-    # discussion with Alan has concluded that it is needed because 
-    # the values are subsequently used in the computation of a weighted average
-    interval_prob <- interval_prob / matrix(colSums(interval_prob), nrow = kn - 1, ncol = ncol(interval_prob), byrow = T)
-
-    interval_prob[is.na(interval_prob)] <- 0
+    # Convert the column of vectors back to a matrix
+    interval_prob <- matrix(unlist(interval_prob_vec), nrow = kn - 1, ncol = length(SODFreq), byrow = F)
 
 
     #######################
@@ -251,6 +208,8 @@ PArisk <- function(
     # Total annual time spent intoxicated over all levels of consumption
     Time_intox_sum <- colSums(Time_intox)
 
+    rm(Duration_m, interval_prob, Duration_h, x, y)
+
   }
 
   #######################
@@ -355,14 +314,8 @@ PArisk <- function(
       FUN.VALUE = numeric(1)
     )
 
-    # 
-    # rm(
-    #   grams_ethanol,
-    #   v, v1, logitp, p
-    # )
-    # gc()
-    # 
-    # 
+
+
 
 
     return(Annual_risk)