From 3127ae3e19514f4a6a6b67ba25fc28d1141e7c71 Mon Sep 17 00:00:00 2001 From: Davide Garolini Date: Tue, 12 Sep 2023 16:35:25 +0200 Subject: [PATCH] Change to poisson dist to see if it is the glm fit to be unstable (#1059) This is not a fix. We need tomorrow's integration tests to see if also the fit changes in the integration tests --- tests/testthat/_snaps/summarize_glm_count.md | 80 +++++++++++++++++-- ... => utils_default_stats_formats_labels.md} | 0 tests/testthat/test-summarize_glm_count.R | 22 +++-- 3 files changed, 90 insertions(+), 12 deletions(-) rename tests/testthat/_snaps/{utils_defaults_handling.md => utils_default_stats_formats_labels.md} (100%) diff --git a/tests/testthat/_snaps/summarize_glm_count.md b/tests/testthat/_snaps/summarize_glm_count.md index a2e15c0bd7..eeef2e6e7c 100644 --- a/tests/testthat/_snaps/summarize_glm_count.md +++ b/tests/testthat/_snaps/summarize_glm_count.md @@ -65,13 +65,79 @@ # h_ppmeans works with healthy input + Code + fits + Output + $glm_fit + + Call: stats::glm(formula = formula, family = stats::poisson(link = "log"), + data = .df_row, offset = offset) + + Coefficients: + (Intercept) REGION1Asia REGION1Eurasia + 2.01066 0.07631 0.64426 + REGION1Europe REGION1North America REGION1South America + 2.13097 -0.07450 0.38102 + ARMCDARM B ARMCDARM C + 0.11048 -0.17694 + + Degrees of Freedom: 199 Total (i.e. Null); 192 Residual + Null Deviance: 983.8 + Residual Deviance: 939 AIC: 1498 + + $emmeans_fit + ARMCD rate SE df asymp.LCL asymp.UCL + ARM A 12.6 1.238 Inf 10.43 15.3 + ARM B 14.1 1.285 Inf 11.81 16.9 + ARM C 10.6 0.971 Inf 8.85 12.7 + + Results are averaged over the levels of: REGION1 + Confidence level used: 0.95 + Intervals are back-transformed from the log scale + + +--- + + Code + fits2 + Output + $glm_fit + + Call: stats::glm(formula = formula, family = stats::quasipoisson(link = "log"), + data = .df_row, offset = offset) + + Coefficients: + (Intercept) REGION1Asia REGION1Eurasia + 2.01066 0.07631 0.64426 + REGION1Europe REGION1North America REGION1South America + 2.13097 -0.07450 0.38102 + ARMCDARM B ARMCDARM C + 0.11048 -0.17694 + + Degrees of Freedom: 199 Total (i.e. Null); 192 Residual + Null Deviance: 983.8 + Residual Deviance: 939 AIC: NA + + $emmeans_fit + ARMCD rate SE df asymp.LCL asymp.UCL + ARM A 12.6 5.20 Inf 5.65 28.3 + ARM B 14.1 5.39 Inf 6.68 29.8 + ARM C 10.6 4.07 Inf 4.98 22.5 + + Results are averaged over the levels of: REGION1 + Confidence level used: 0.95 + Intervals are back-transformed from the log scale + + +--- + Code res Output rate asymp.LCL asymp.UCL ARM - A: Drug X 3.07 2.202774 4.278651 A: Drug X - B: Placebo 3.07 2.202774 4.278651 B: Placebo - C: Combination 3.07 2.202774 4.278651 C: Combination + A: Drug X 3.07 2.836527 3.32269 A: Drug X + B: Placebo 3.07 2.836527 3.32269 B: Placebo + C: Combination 3.07 2.836527 3.32269 C: Combination # s_glm_count works with healthy input @@ -87,7 +153,7 @@ [1] "Adjusted Rate" $rate_ci - [1] 1.983340 6.127155 + [1] 3.047667 3.987387 attr(,"label") [1] "95% CI" @@ -121,7 +187,7 @@ [1] "Adjusted Rate" $rate_ci - [1] 1.983340 6.127155 + [1] 3.047667 3.987387 attr(,"label") [1] "95% CI" @@ -131,12 +197,12 @@ [1] "Adjusted Rate Ratio" $rate_ratio_ci - [1] 0.3974979 0.3067371 2.0169939 1.8347687 + [1] 0.7378778 0.6062152 1.0865633 0.9283695 attr(,"label") [1] "95% CI" $pval - [1] 0.7897470 0.5287652 + [1] 0.263119218 0.008203621 attr(,"label") [1] "p-value" diff --git a/tests/testthat/_snaps/utils_defaults_handling.md b/tests/testthat/_snaps/utils_default_stats_formats_labels.md similarity index 100% rename from tests/testthat/_snaps/utils_defaults_handling.md rename to tests/testthat/_snaps/utils_default_stats_formats_labels.md diff --git a/tests/testthat/test-summarize_glm_count.R b/tests/testthat/test-summarize_glm_count.R index 4ced0c8070..9d7eb2d94d 100644 --- a/tests/testthat/test-summarize_glm_count.R +++ b/tests/testthat/test-summarize_glm_count.R @@ -169,11 +169,23 @@ testthat::test_that("h_ppmeans works with healthy input", { anl$AVAL_f <- as.factor(anl$AVAL) fits <- h_glm_count( + .var = "AVAL", + .df_row = anl, + variables = list(arm = "ARMCD", offset = "lgTMATRSK", covariates = c("REGION1")), + distribution = "poisson" + ) + + testthat::expect_snapshot(fits) + + fits2 <- h_glm_count( .var = "AVAL", .df_row = anl, variables = list(arm = "ARMCD", offset = "lgTMATRSK", covariates = c("REGION1")), distribution = "quasipoisson" ) + + testthat::expect_snapshot(fits2) + result <- h_ppmeans( obj = fits$glm_fit, .df_row = anl, @@ -182,7 +194,7 @@ testthat::test_that("h_ppmeans works with healthy input", { ) res <- testthat::expect_silent(result) - testthat::expect_snapshot(res) + testthat::expect_snapshot(res) # diff }) testthat::test_that("s_glm_count works with healthy input", { @@ -198,12 +210,12 @@ testthat::test_that("s_glm_count works with healthy input", { .in_ref_col = TRUE, variables = list(arm = "ARMCD", offset = "lgTMATRSK", covariates = c("REGION1")), conf_level = 0.95, - distribution = "quasipoisson", + distribution = "poisson", rate_mean_method = "ppmeans" ) res <- testthat::expect_silent(result) - testthat::expect_snapshot(res) + testthat::expect_snapshot(res) # diff }) testthat::test_that("s_glm_count works with no reference group selected.", { @@ -221,12 +233,12 @@ testthat::test_that("s_glm_count works with no reference group selected.", { filter(ARMCD == "ARM B"), variables = list(arm = "ARMCD", offset = "lgTMATRSK", covariates = c("REGION1")), conf_level = 0.95, - distribution = "quasipoisson", + distribution = "poisson", rate_mean_method = "ppmeans" ) res <- testthat::expect_silent(result) - testthat::expect_snapshot(res) + testthat::expect_snapshot(res) # diff }) testthat::test_that("s_glm_count fails wrong inputs", {