fixed bug in generate_data and adapted tests (#29)

DESCRIPTION

@@ -1,8 +1,8 @@
 Type: Package
 Package: savvyr
 Title: Survival Analysis for AdVerse Events with VarYing Follow-Up Times
-Version: 0.1.0
-Date: 2024-02-20
+Version: 0.1.1
+Date: 2024-06-06
 Authors@R: c(
     person("Thomas", "Kuenzel", email = "thomas.kuenzel@roche.com", role = c("aut", "cre")),
     person("Kaspar", "Rufibach", email = "kaspar.rufibach@roche.com", role = "aut"),

NEWS.md

@@ -1,3 +1,9 @@
+# savvyr 0.1.1
+
+### Bug Fixes
+
+- Changed the way data is generated in function 'generate_data': The variable 'type_event' is now generated such that it accurately reflects the specified hazards for AE, death and soft competing events.
+
 # savvyr 0.1.0
 
 - First CRAN version of the package.
@@ -7,3 +13,6 @@
 
 - Estimators that do not account for competing events (incidence proportion, incidence density, Inverse Kaplan Meier).
 - Estimators accounting for competing events (incidence proportion accounting for competing events and Aalen-Johansen, both first with death only as hard competing event, or using all competing events).
+
+
+

R/generate_data.R

@@ -39,7 +39,15 @@ generate_data <- function(n,
                           haz_death,
                           haz_soft) {
   assert_count(n, positive = TRUE)
-  assert_numeric(cens, lower = 0, finite = TRUE, any.missing = FALSE, len = 2L, unique = TRUE, sorted = TRUE)
+  assert_numeric(
+    cens,
+    lower = 0,
+    finite = TRUE,
+    any.missing = FALSE,
+    len = 2L,
+    unique = TRUE,
+    sorted = TRUE
+  )
   assert_number(haz_ae, finite = TRUE)
   assert_number(haz_death, finite = TRUE)
   assert_number(haz_soft, finite = TRUE)
@@ -52,11 +60,14 @@ generate_data <- function(n,
   )
   haz_all <- sum(haz)
   result$time_to_event <- stats::rexp(n = n, rate = haz_all)
-  result$type_of_event <- 1L + stats::rbinom(
-    n = n,
-    size = 2,
-    prob = haz / haz_all
+
+  result$type_of_event <- sample(
+    1:3,
+    size = n,
+    prob = haz / haz_all,
+    replace = TRUE
   )
+
   result$cens <- stats::runif(n = n, min = cens[1L], max = cens[2L])
   result$type_of_event <- ifelse(
     result$time_to_event <= result$cens,

tests/testthat/test-aalen_johansen.R

@@ -8,7 +8,7 @@ test_that("Aalen Johansen works as expected", {
     haz_soft = 0.5
   )
   result <- aalen_johansen(data = df, ce = 2, tau = 4)
-  expected <- c(ae_prob = 0.2719, ae_prob_var = 0.0119, ce_prob = 0.7281, ce_prob_var = 0.0119)
+  expected <- c(ae_prob = 0.4792, ae_prob_var = 0.0266, ce_prob = 0.5208, ce_prob_var = 0.0266)
   expect_equal(result, expected, tolerance = 1e-4)
 })
 
@@ -39,7 +39,7 @@ test_that("Aalen Johansen works without competing events", {
   df <- df[df$type_of_event != 2, ]
   df <- df[df$type_of_event != 3, ]
   result <- aalen_johansen(data = df, ce = 2, tau = 4)
-  expected <- c(ae_prob = 0.5897, ae_prob_var = 0.0404, ce_prob = 0, ce_prob_var = 0)
+  expected <- c(ae_prob = 0.7643, ae_prob_var = 0.0363, ce_prob = 0, ce_prob_var = 0)
   expect_equal(result, expected, tolerance = 1e-4)
 })
 

tests/testthat/test-generate_data.R

@@ -12,3 +12,17 @@ test_that("generate_data works as expected", {
   expect_identical(result$id, 1:10)
   expect_integer(result$type_of_event, lower = 0, upper = 3, any.missing = FALSE)
 })
+
+test_that("generate_data correctly ... ", {
+  set.seed(123)
+  df <- generate_data(
+    n = 10^6,
+    cens = c(1000, 1001),
+    haz_ae = 0.2,
+    haz_death = 0.3,
+    haz_soft = 0.5
+  )
+  result <- inc_prop(data = df, tau = 4)
+  expected <- c(ae_prob = 0.2)
+  expect_equal(result["ae_prob"], expected, tolerance = 0.1)
+})

tests/testthat/test-inc_prop.R

@@ -8,6 +8,6 @@ test_that("inc_prop works as expected", {
     haz_soft = 0.5
   )
   result <- inc_prop(data = df, tau = 4)
-  expected <- c(ae_prob = 0.2, ae_prob_var = 0.0064)
+  expected <- c(ae_prob = 0.24, ae_prob_var = 0.0073)
   expect_equal(result, expected, tolerance = 1e-4)
 })

tests/testthat/test-one_minus_kaplan_meier.R

@@ -8,7 +8,7 @@ test_that("one_minus_kaplan_meier works as expected", {
     haz_soft = 0.5
   )
   result <- one_minus_kaplan_meier(data = df, tau = 4)
-  expected <- c(ae_prob = 0.3771350, ae_prob_var = 0.0260535)
+  expected <- c(ae_prob = 0.5865, ae_prob_var = 0.0394)
   expect_equal(result, expected, tolerance = 1e-4)
 })
 

tests/testthat/test-prop_trans_inc_dens.R

@@ -8,6 +8,6 @@ test_that("prop_trans_inc_dens works as expected", {
     haz_soft = 1.2
   )
   result <- prop_trans_inc_dens(data = df, tau = 0.1)
-  expected <- c(ae_prob = 0.33325, ae_prob_var = 0.0365)
+  expected <- c(ae_prob = 0.5554, ae_prob_var = 0.0325)
   expect_equal(result, expected, tolerance = 1e-4)
 })

tests/testthat/test-prop_trans_inc_dens_ce.R

@@ -8,6 +8,6 @@ test_that("prop_trans_inc_dens_ce works as expected", {
     haz_soft = 1.2
   )
   result <- prop_trans_inc_dens_ce(data = df, ce = 2, tau = 0.1)
-  expected <- c(ae_prob = 0.296148137, ae_prob_var = 0.002755189)
+  expected <- c(ae_prob = 0.18114, ae_prob_var = 0.00205)
   expect_equal(result, expected, tolerance = 1e-4)
 })