From b4ba7467ad67350db6deba5f609179e6fa0f3d28 Mon Sep 17 00:00:00 2001 From: Daniel Date: Sun, 24 Sep 2023 18:22:53 +0200 Subject: [PATCH 01/14] close #524 --- R/check_outliers.R | 19 +++++++++++++++++-- 1 file changed, 17 insertions(+), 2 deletions(-) diff --git a/R/check_outliers.R b/R/check_outliers.R index c1af3bb1f..519ab64e0 100644 --- a/R/check_outliers.R +++ b/R/check_outliers.R @@ -167,7 +167,9 @@ #' value for outliers classification. Refer to the help-file of #' [ICSOutlier::ics.outlier()] to get more details about this procedure. #' Note that `method = "ics"` requires both **ICS** and **ICSOutlier** -#' to be installed, and that it takes some time to compute the results. +#' to be installed, and that it takes some time to compute the results. You +#' can speed up computation time using parallel computing. Set the number of +#' cores to use with `options(mc.cores = 4)` (for example). #' #' - **OPTICS**: #' The Ordering Points To Identify the Clustering Structure (OPTICS) algorithm @@ -1764,7 +1766,20 @@ check_outliers.metabin <- check_outliers.metagen n_cores <- if (!requireNamespace("parallel", quietly = TRUE)) { NULL } else { - max(1L, parallel::detectCores() - 2L, na.rm = TRUE) + getOption("mc.cores", 1L) + } + + # tell user about n-cores option + if (is.null(n_cores)) { + insight::format_alert( + "Package `parallel` is not installed. `check_outliers()` will run on a single core.", + "Install package `parallel` and set, for example, `options(mc.cores = 4)` to run on multiple cores." + ) + } else if (n_cores == 1) { + insight::format_alert( + "Package `parallel` is installed, but `check_outliers()` will run on a single core.", + "To use multiple cores, set `options(mc.cores = 4)` (for example)." + ) } # Run algorithm From 692860a5702b47245d02bb051a170a5651462f24 Mon Sep 17 00:00:00 2001 From: Daniel Date: Sun, 24 Sep 2023 18:23:22 +0200 Subject: [PATCH 02/14] rd --- man/check_outliers.Rd | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/man/check_outliers.Rd b/man/check_outliers.Rd index e3c218b7a..bb95fadf8 100644 --- a/man/check_outliers.Rd +++ b/man/check_outliers.Rd @@ -195,7 +195,9 @@ of 0.025 (corresponding to the 2.5\\% most extreme observations) as a cut-off value for outliers classification. Refer to the help-file of \code{\link[ICSOutlier:ics.outlier]{ICSOutlier::ics.outlier()}} to get more details about this procedure. Note that \code{method = "ics"} requires both \strong{ICS} and \strong{ICSOutlier} -to be installed, and that it takes some time to compute the results. +to be installed, and that it takes some time to compute the results. You +can speed up computation time using parallel computing. Set the number of +cores to use with \code{options(mc.cores = 4)} (for example). \item \strong{OPTICS}: The Ordering Points To Identify the Clustering Structure (OPTICS) algorithm (Ankerst et al., 1999) is using similar concepts to DBSCAN (an unsupervised From 7f4719ab8f9662599c495157f4dede249a84dcd2 Mon Sep 17 00:00:00 2001 From: Daniel Date: Sun, 24 Sep 2023 18:25:07 +0200 Subject: [PATCH 03/14] news, dev --- DESCRIPTION | 2 +- NEWS.md | 7 +++++++ 2 files changed, 8 insertions(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index fa29f37d4..4853a45d3 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Type: Package Package: performance Title: Assessment of Regression Models Performance -Version: 0.10.5.2 +Version: 0.10.5.3 Authors@R: c(person(given = "Daniel", family = "Lüdecke", diff --git a/NEWS.md b/NEWS.md index 857924e29..40e57d346 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,12 @@ # performance (development version) +## Changes to functions + +* `check_outliers()` for method `"ics"` now detects number of available cores + for parallel computing via the `"mc.cores"` option. This is more robust than + the previous method, which used `parallel::detectCores()`. Now you should + set the number of cores via `options(mc.cores = 4)`. + # performance 0.10.5 ## Changes to functions From 462e8aef4757f628a154c2752ac914b2f7d1be06 Mon Sep 17 00:00:00 2001 From: Daniel Date: Sun, 24 Sep 2023 18:59:59 +0200 Subject: [PATCH 04/14] docs --- R/check_outliers.R | 8 +++----- man/check_outliers.Rd | 7 +++---- 2 files changed, 6 insertions(+), 9 deletions(-) diff --git a/R/check_outliers.R b/R/check_outliers.R index 519ab64e0..e98e037c0 100644 --- a/R/check_outliers.R +++ b/R/check_outliers.R @@ -301,6 +301,7 @@ #' group_iris <- datawizard::data_group(iris, "Species") #' check_outliers(group_iris) #' +#' @examplesIf require("see") && require("bigutilsr") && require("loo") && require("MASS") && require("ICSOutlier") && require("ICS") && require("dbscan") #' \donttest{ #' # You can also run all the methods #' check_outliers(data, method = "all") @@ -317,10 +318,7 @@ #' model <- lm(disp ~ mpg + hp, data = mt2) #' #' outliers_list <- check_outliers(model) -#' -#' if (require("see")) { -#' plot(outliers_list) -#' } +#' plot(outliers_list) #' #' insight::get_data(model)[outliers_list, ] # Show outliers data #' } @@ -508,7 +506,7 @@ check_outliers.default <- function(x, num.df <- outlier_count$all[!names(outlier_count$all) %in% c("Row", ID)] if (isTRUE(nrow(num.df) > 0)) { - num.df <- datawizard::change_code( + num.df <- datawizard::recode_values( num.df, recode = list(`2` = "(Multivariate)") ) diff --git a/man/check_outliers.Rd b/man/check_outliers.Rd index bb95fadf8..c19d4ecb0 100644 --- a/man/check_outliers.Rd +++ b/man/check_outliers.Rd @@ -292,6 +292,7 @@ filtered_data <- data[outliers_info$Outlier < 0.1, ] group_iris <- datawizard::data_group(iris, "Species") check_outliers(group_iris) +\dontshow{if (require("see") && require("bigutilsr") && require("loo") && require("MASS") && require("ICSOutlier") && require("ICS") && require("dbscan")) (if (getRversion() >= "3.4") withAutoprint else force)(\{ # examplesIf} \donttest{ # You can also run all the methods check_outliers(data, method = "all") @@ -308,13 +309,11 @@ mt2 <- rbind(mt1, data.frame( model <- lm(disp ~ mpg + hp, data = mt2) outliers_list <- check_outliers(model) - -if (require("see")) { - plot(outliers_list) -} +plot(outliers_list) insight::get_data(model)[outliers_list, ] # Show outliers data } +\dontshow{\}) # examplesIf} } \references{ \itemize{ From ddaebc76bce529c9b926bcc144c0cc045c4dd928 Mon Sep 17 00:00:00 2001 From: Daniel Date: Wed, 27 Sep 2023 15:08:20 +0200 Subject: [PATCH 05/14] Recursion in `logLik.plm()` that slows down computations: (#623) Co-authored-by: Indrajeet Patil --- R/logLik.R | 2 +- tests/testthat/test-logLik.R | 27 +++++++++++++++++++++++++++ 2 files changed, 28 insertions(+), 1 deletion(-) create mode 100644 tests/testthat/test-logLik.R diff --git a/R/logLik.R b/R/logLik.R index 90c12930b..5c42f41a5 100644 --- a/R/logLik.R +++ b/R/logLik.R @@ -55,7 +55,7 @@ logLik.plm <- function(object, ...) { attr(val, "nall") <- N0 attr(val, "nobs") <- N - attr(val, "df") <- insight::get_df(object, type = "model") + attr(val, "df") <- insight::n_parameters(object) + 1L class(val) <- "logLik" val diff --git a/tests/testthat/test-logLik.R b/tests/testthat/test-logLik.R new file mode 100644 index 000000000..2691f1d3c --- /dev/null +++ b/tests/testthat/test-logLik.R @@ -0,0 +1,27 @@ +test_that("logLik", { + skip_if_not_installed("plm") + skip_if_not_installed("withr") + + withr::local_options(list(expressions = 25)) + set.seed(1) + nnn <- 100 + ddta <- data.frame( + x1 = rnorm(nnn), + x2 = rnorm(nnn), + id = rep_len(1:(nnn / 10), nnn), + year = rep_len(1:11, nnn) + ) + ddta$y <- ddta$x1 * 0.5 - ddta$x2 * 0.5 + rnorm(nnn) + + m1 <- lm(y ~ x1 + x2, data = ddta) + l1 <- logLik(m1) + + m2 <- plm( + y ~ x1 + x2, + data = ddta, + model = "pooling", + index = c("id", "year") + ) + l2 <- logLik(m2) + expect_equal(l1, l2, tolerance = 1e-3, ignore_attr = TRUE) +}) From aca7e9964cba8a6625fe22e0c3b3a6eb00308e5a Mon Sep 17 00:00:00 2001 From: Indrajeet Patil Date: Wed, 27 Sep 2023 19:32:24 +0530 Subject: [PATCH 06/14] Experiment with autobump deps PR workflow (#624) --- .../workflows/update-to-latest-easystats.yml | 54 +++++++++++++++++++ 1 file changed, 54 insertions(+) create mode 100644 .github/workflows/update-to-latest-easystats.yml diff --git a/.github/workflows/update-to-latest-easystats.yml b/.github/workflows/update-to-latest-easystats.yml new file mode 100644 index 000000000..90f6f3fc3 --- /dev/null +++ b/.github/workflows/update-to-latest-easystats.yml @@ -0,0 +1,54 @@ +name: update-to-latest-easystats + +on: + push: + branches: [main, master] + pull_request: + branches: [main, master] + +# Run on a schedule (once a month) +# schedule: +# - cron: "0 0 1 * *" + +jobs: + update-to-latest-easystats: + runs-on: ubuntu-latest + steps: + - name: Checkout code + uses: actions/checkout@v4 + + - name: Set up R + uses: r-lib/actions/setup-r@v2 + + - uses: r-lib/actions/setup-r-dependencies@v2 + with: + dependencies: '"hard"' + extra-packages: | + any::usethis + + - name: Update DESCRIPTION file + id: update_description + run: | + usethis::use_latest_dependencies() + shell: Rscript {0} + + - name: Check for changes in DESCRIPTION + id: check_description + run: | + git status + git add --all + + - name: Create pull request + #if: ${{ steps.update_description.outputs.changed }} + uses: peter-evans/create-pull-request@v5 + with: + token: ${{ secrets.GITHUB_TOKEN }} + base: main + branch: desc-${{ github.ref_name }}-${{ github.job }} + branch-suffix: timestamp + delete-branch: true + title: "Update `DESCRIPTION` to use latest 'easystats' dependencies" + body: "Automatically updated the `DESCRIPTION` file using `usethis::use_latest_dependencies()`." + labels: "auto-update" + add-paths: | + DESCRIPTION From 4870b8c166338a6424b90acd34b0f2b85c486af6 Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" <41898282+github-actions[bot]@users.noreply.github.com> Date: Wed, 27 Sep 2023 19:32:36 +0530 Subject: [PATCH 07/14] Update `DESCRIPTION` to use latest 'easystats' dependencies (#625) Co-authored-by: IndrajeetPatil --- DESCRIPTION | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 4853a45d3..302314eb1 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -69,9 +69,9 @@ BugReports: https://github.com/easystats/performance/issues Depends: R (>= 3.6) Imports: - bayestestR (>= 0.13.0), - insight (>= 0.19.4), - datawizard (>= 0.7.0), + bayestestR (>= 0.13.1), + insight (>= 0.19.5), + datawizard (>= 0.9.0), methods, stats, utils From a5647f900f91773f387353ac4d856b3f0221c63b Mon Sep 17 00:00:00 2001 From: Indrajeet Patil Date: Wed, 27 Sep 2023 16:20:55 +0200 Subject: [PATCH 08/14] bump version --- DESCRIPTION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index 302314eb1..3b95351d9 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Type: Package Package: performance Title: Assessment of Regression Models Performance -Version: 0.10.5.3 +Version: 0.10.5.4 Authors@R: c(person(given = "Daniel", family = "Lüdecke", From 3be01281a236ca945a71118d14dfa026739c25be Mon Sep 17 00:00:00 2001 From: Indrajeet Patil Date: Wed, 27 Sep 2023 16:31:37 +0200 Subject: [PATCH 09/14] expt over [skip ci] --- .../workflows/update-to-latest-easystats.yml | 54 ------------------- 1 file changed, 54 deletions(-) delete mode 100644 .github/workflows/update-to-latest-easystats.yml diff --git a/.github/workflows/update-to-latest-easystats.yml b/.github/workflows/update-to-latest-easystats.yml deleted file mode 100644 index 90f6f3fc3..000000000 --- a/.github/workflows/update-to-latest-easystats.yml +++ /dev/null @@ -1,54 +0,0 @@ -name: update-to-latest-easystats - -on: - push: - branches: [main, master] - pull_request: - branches: [main, master] - -# Run on a schedule (once a month) -# schedule: -# - cron: "0 0 1 * *" - -jobs: - update-to-latest-easystats: - runs-on: ubuntu-latest - steps: - - name: Checkout code - uses: actions/checkout@v4 - - - name: Set up R - uses: r-lib/actions/setup-r@v2 - - - uses: r-lib/actions/setup-r-dependencies@v2 - with: - dependencies: '"hard"' - extra-packages: | - any::usethis - - - name: Update DESCRIPTION file - id: update_description - run: | - usethis::use_latest_dependencies() - shell: Rscript {0} - - - name: Check for changes in DESCRIPTION - id: check_description - run: | - git status - git add --all - - - name: Create pull request - #if: ${{ steps.update_description.outputs.changed }} - uses: peter-evans/create-pull-request@v5 - with: - token: ${{ secrets.GITHUB_TOKEN }} - base: main - branch: desc-${{ github.ref_name }}-${{ github.job }} - branch-suffix: timestamp - delete-branch: true - title: "Update `DESCRIPTION` to use latest 'easystats' dependencies" - body: "Automatically updated the `DESCRIPTION` file using `usethis::use_latest_dependencies()`." - labels: "auto-update" - add-paths: | - DESCRIPTION From 45cf8d0861ec0f9b84d18c4a5eecce66b6e0818d Mon Sep 17 00:00:00 2001 From: Indrajeet Patil Date: Wed, 27 Sep 2023 16:54:56 +0200 Subject: [PATCH 10/14] Create update-to-latest-easystats.yaml --- .github/workflows/update-to-latest-easystats.yaml | 10 ++++++++++ 1 file changed, 10 insertions(+) create mode 100644 .github/workflows/update-to-latest-easystats.yaml diff --git a/.github/workflows/update-to-latest-easystats.yaml b/.github/workflows/update-to-latest-easystats.yaml new file mode 100644 index 000000000..7718d763a --- /dev/null +++ b/.github/workflows/update-to-latest-easystats.yaml @@ -0,0 +1,10 @@ +on: + schedule: + # Check for dependency updates once a month + - cron: "0 0 1 * *" + +name: update-to-latest-easystats + +jobs: + update-to-latest-easystats: + uses: easystats/workflows/.github/workflows/update-to-latest-easystats.yaml@main From 39e098b8cc915295e7d5c51a473fbd8b4e44249b Mon Sep 17 00:00:00 2001 From: Indrajeet Patil Date: Thu, 28 Sep 2023 13:10:51 +0200 Subject: [PATCH 11/14] Update html-5-check.yaml --- .github/workflows/html-5-check.yaml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/.github/workflows/html-5-check.yaml b/.github/workflows/html-5-check.yaml index 1439a3228..f25b4eeaf 100644 --- a/.github/workflows/html-5-check.yaml +++ b/.github/workflows/html-5-check.yaml @@ -6,8 +6,8 @@ on: pull_request: branches: [main, master] -name: HTML5 check +name: html-5-check jobs: - HTML5-check: + html-5-check: uses: easystats/workflows/.github/workflows/html-5-check.yaml@main From 05336ec12f75b5eb696f698c97028ae46dee0b00 Mon Sep 17 00:00:00 2001 From: Daniel Date: Fri, 29 Sep 2023 13:30:12 +0200 Subject: [PATCH 12/14] Error: check_model() not implemented for models of class lm yet (#630) --- DESCRIPTION | 2 +- NEWS.md | 5 +++++ R/check_model_diagnostics.R | 8 ++++---- 3 files changed, 10 insertions(+), 5 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 3b95351d9..19bada6a3 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Type: Package Package: performance Title: Assessment of Regression Models Performance -Version: 0.10.5.4 +Version: 0.10.5.5 Authors@R: c(person(given = "Daniel", family = "Lüdecke", diff --git a/NEWS.md b/NEWS.md index 40e57d346..47fad4e59 100644 --- a/NEWS.md +++ b/NEWS.md @@ -7,6 +7,11 @@ the previous method, which used `parallel::detectCores()`. Now you should set the number of cores via `options(mc.cores = 4)`. +## Bug fixes + +* Fixed issues is `check_model()` for models that used data sets with + variables of class `"haven_labelled"`. + # performance 0.10.5 ## Changes to functions diff --git a/R/check_model_diagnostics.R b/R/check_model_diagnostics.R index 3b86799cc..371e32609 100644 --- a/R/check_model_diagnostics.R +++ b/R/check_model_diagnostics.R @@ -153,7 +153,7 @@ # prepare data for normality of residuals plot ---------------------------------- .diag_norm <- function(model, verbose = TRUE) { - r <- try(stats::residuals(model), silent = TRUE) + r <- try(as.numeric(stats::residuals(model)), silent = TRUE) if (inherits(r, "try-error")) { insight::format_alert(sprintf("Non-normality of residuals could not be computed. Cannot extract residuals from objects of class '%s'.", class(model)[1])) @@ -183,7 +183,7 @@ n_params <- tryCatch(model$rank, error = function(e) insight::n_parameters(model)) infl <- stats::influence(model, do.coef = FALSE) - resid <- insight::get_residuals(model) + resid <- as.numeric(insight::get_residuals(model)) std_resid <- tryCatch(stats::rstandard(model, infl), error = function(e) resid) @@ -212,8 +212,8 @@ ncv <- tryCatch( { data.frame( - x = stats::fitted(model), - y = stats::residuals(model) + x = as.numeric(stats::fitted(model)), + y = as.numeric(stats::residuals(model)) ) }, error = function(e) { From 198578d01e4e85cb1907761b861a59f57ef945bd Mon Sep 17 00:00:00 2001 From: Daniel Date: Fri, 29 Sep 2023 17:34:27 +0200 Subject: [PATCH 13/14] lintr --- R/check_model_diagnostics.R | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/R/check_model_diagnostics.R b/R/check_model_diagnostics.R index 371e32609..46a64f9c2 100644 --- a/R/check_model_diagnostics.R +++ b/R/check_model_diagnostics.R @@ -174,10 +174,10 @@ if (inherits(model, "lm", which = TRUE) == 1) { cook_levels <- round(stats::qf(0.5, s$fstatistic[2], s$fstatistic[3]), 2) - } else if (!is.null(threshold)) { - cook_levels <- threshold - } else { + } else if (is.null(threshold)) { cook_levels <- c(0.5, 1) + } else { + cook_levels <- threshold } n_params <- tryCatch(model$rank, error = function(e) insight::n_parameters(model)) From 32199f27f0caa3016ba0a60d4bc1eeff40c339e6 Mon Sep 17 00:00:00 2001 From: Daniel Date: Wed, 4 Oct 2023 14:10:13 +0200 Subject: [PATCH 14/14] Support for nestedLogit models (#633) --- DESCRIPTION | 5 +- NAMESPACE | 6 ++ NEWS.md | 4 + R/check_itemscale.R | 8 +- R/model_performance.lm.R | 22 ++++- R/r2.R | 81 +++++++++---------- R/r2_coxsnell.R | 20 +++++ R/r2_nagelkerke.R | 20 +++++ R/r2_tjur.R | 31 +++++++ man/check_itemscale.Rd | 2 - tests/testthat/_snaps/check_distribution.md | 21 +++++ tests/testthat/_snaps/nestedLogit.md | 12 +++ tests/testthat/test-binned_residuals.R | 67 +++++++++++++++ tests/testthat/test-check_autocorrelation.R | 7 ++ tests/testthat/test-check_distribution.R | 36 +++++++++ .../testthat/test-check_heterogeneity_bias.R | 33 ++++++++ tests/testthat/test-nestedLogit.R | 66 +++++++++++++++ 17 files changed, 387 insertions(+), 54 deletions(-) create mode 100644 tests/testthat/_snaps/check_distribution.md create mode 100644 tests/testthat/_snaps/nestedLogit.md create mode 100644 tests/testthat/test-binned_residuals.R create mode 100644 tests/testthat/test-check_autocorrelation.R create mode 100644 tests/testthat/test-check_distribution.R create mode 100644 tests/testthat/test-check_heterogeneity_bias.R create mode 100644 tests/testthat/test-nestedLogit.R diff --git a/DESCRIPTION b/DESCRIPTION index 19bada6a3..933df479a 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Type: Package Package: performance Title: Assessment of Regression Models Performance -Version: 0.10.5.5 +Version: 0.10.5.6 Authors@R: c(person(given = "Daniel", family = "Lüdecke", @@ -117,6 +117,7 @@ Suggests: mgcv, mlogit, multimode, + nestedLogit, nlme, nonnest2, ordinal, @@ -149,4 +150,4 @@ Config/Needs/website: r-lib/pkgdown, easystats/easystatstemplate Config/rcmdcheck/ignore-inconsequential-notes: true -Remotes: easystats/see, easystats/parameters +Remotes: easystats/see, easystats/parameters, easystats/insight diff --git a/NAMESPACE b/NAMESPACE index b5c3ee92b..e573d5ce0 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -180,6 +180,7 @@ S3method(model_performance,model_fit) S3method(model_performance,multinom) S3method(model_performance,negbinirr) S3method(model_performance,negbinmfx) +S3method(model_performance,nestedLogit) S3method(model_performance,plm) S3method(model_performance,poissonirr) S3method(model_performance,poissonmfx) @@ -371,6 +372,7 @@ S3method(r2,model_fit) S3method(r2,multinom) S3method(r2,negbinirr) S3method(r2,negbinmfx) +S3method(r2,nestedLogit) S3method(r2,ols) S3method(r2,phylolm) S3method(r2,plm) @@ -413,6 +415,7 @@ S3method(r2_coxsnell,mclogit) S3method(r2_coxsnell,multinom) S3method(r2_coxsnell,negbinirr) S3method(r2_coxsnell,negbinmfx) +S3method(r2_coxsnell,nestedLogit) S3method(r2_coxsnell,poissonirr) S3method(r2_coxsnell,poissonmfx) S3method(r2_coxsnell,polr) @@ -468,6 +471,7 @@ S3method(r2_nagelkerke,mclogit) S3method(r2_nagelkerke,multinom) S3method(r2_nagelkerke,negbinirr) S3method(r2_nagelkerke,negbinmfx) +S3method(r2_nagelkerke,nestedLogit) S3method(r2_nagelkerke,poissonirr) S3method(r2_nagelkerke,poissonmfx) S3method(r2_nagelkerke,polr) @@ -479,6 +483,8 @@ S3method(r2_posterior,BFBayesFactor) S3method(r2_posterior,brmsfit) S3method(r2_posterior,stanmvreg) S3method(r2_posterior,stanreg) +S3method(r2_tjur,default) +S3method(r2_tjur,nestedLogit) S3method(residuals,BFBayesFactor) S3method(residuals,check_normality_numeric) S3method(residuals,iv_robust) diff --git a/NEWS.md b/NEWS.md index 47fad4e59..06eee8f47 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,9 @@ # performance (development version) +## General + +* Support for `nestedLogit` models. + ## Changes to functions * `check_outliers()` for method `"ics"` now detects number of available cores diff --git a/R/check_itemscale.R b/R/check_itemscale.R index 8d8b71082..5aa704618 100644 --- a/R/check_itemscale.R +++ b/R/check_itemscale.R @@ -40,8 +40,6 @@ #' - Briggs SR, Cheek JM (1986) The role of factor analysis in the development #' and evaluation of personality scales. Journal of Personality, 54(1), #' 106-148. doi: 10.1111/j.1467-6494.1986.tb00391.x -#' - Trochim WMK (2008) Types of Reliability. -#' ([web](https://conjointly.com/kb/types-of-reliability/)) #' #' @examplesIf require("parameters") && require("psych") #' # data generation from '?prcomp', slightly modified @@ -82,9 +80,9 @@ check_itemscale <- function(x) { Mean = vapply(items, mean, numeric(1), na.rm = TRUE), SD = vapply(items, stats::sd, numeric(1), na.rm = TRUE), Skewness = vapply(items, function(i) as.numeric(datawizard::skewness(i)), numeric(1)), - "Difficulty" = item_difficulty(items)$Difficulty, - "Discrimination" = .item_discr, - "alpha if deleted" = .item_alpha, + Difficulty = item_difficulty(items)$Difficulty, + Discrimination = .item_discr, + `alpha if deleted` = .item_alpha, stringsAsFactors = FALSE, check.names = FALSE ) diff --git a/R/model_performance.lm.R b/R/model_performance.lm.R index e55253125..e8e5b07cf 100644 --- a/R/model_performance.lm.R +++ b/R/model_performance.lm.R @@ -125,10 +125,10 @@ model_performance.lm <- function(model, metrics = "all", verbose = TRUE, ...) { if (("LOGLOSS" %in% toupper(metrics)) && isTRUE(info$is_binomial)) { out$Log_loss <- .safe({ .logloss <- performance_logloss(model, verbose = verbose) - if (!is.na(.logloss)) { - .logloss - } else { + if (is.na(.logloss)) { NULL + } else { + .logloss } }) } @@ -253,6 +253,22 @@ model_performance.zeroinfl <- model_performance.lm #' @export model_performance.zerotrunc <- model_performance.lm +#' @export +model_performance.nestedLogit <- function(model, metrics = "all", verbose = TRUE, ...) { + mp <- lapply(model$models, model_performance.lm, metrics = metrics, verbose = verbose, ...) + out <- cbind( + data.frame(Response = names(mp), stringsAsFactors = FALSE), + do.call(rbind, mp) + ) + # need to handle R2 separately + if (any(c("ALL", "R2") %in% toupper(metrics))) { + out$R2 <- unlist(r2_tjur(model)) + } + + row.names(out) <- NULL + class(out) <- unique(c("performance_model", class(out))) + out +} diff --git a/R/r2.R b/R/r2.R index 26a16f9e6..c844f22cd 100644 --- a/R/r2.R +++ b/R/r2.R @@ -49,11 +49,8 @@ r2 <- function(model, ...) { } - - # Default models ----------------------------------------------- - #' @rdname r2 #' @export r2.default <- function(model, ci = NULL, verbose = TRUE, ...) { @@ -115,6 +112,8 @@ r2.lm <- function(model, ci = NULL, ...) { #' @export r2.phylolm <- r2.lm +# helper ------------- + .r2_lm <- function(model_summary, ci = NULL) { out <- list( R2 = model_summary$r.squared, @@ -140,7 +139,6 @@ r2.phylolm <- r2.lm } - #' @export r2.summary.lm <- function(model, ci = NULL, ...) { if (!is.null(ci) && !is.na(ci)) { @@ -150,7 +148,6 @@ r2.summary.lm <- function(model, ci = NULL, ...) { } - #' @export r2.systemfit <- function(model, ...) { out <- lapply(summary(model)$eq, function(model_summary) { @@ -173,11 +170,11 @@ r2.systemfit <- function(model, ...) { #' @export r2.lm_robust <- function(model, ...) { out <- list( - "R2" = tryCatch( + R2 = tryCatch( model[["r.squared"]], error = function(e) NULL ), - "R2_adjusted" = tryCatch( + R2_adjusted = tryCatch( model[["adj.r.squared"]], error = function(e) NULL ) @@ -198,8 +195,6 @@ r2.ols <- function(model, ...) { structure(class = "r2_generic", out) } - - #' @export r2.lrm <- r2.ols @@ -207,7 +202,6 @@ r2.lrm <- r2.ols r2.cph <- r2.ols - #' @export r2.mhurdle <- function(model, ...) { resp <- insight::get_response(model, verbose = FALSE) @@ -230,7 +224,6 @@ r2.mhurdle <- function(model, ...) { } - #' @export r2.aov <- function(model, ci = NULL, ...) { if (!is.null(ci) && !is.na(ci)) { @@ -252,7 +245,6 @@ r2.aov <- function(model, ci = NULL, ...) { } - #' @export r2.mlm <- function(model, ...) { model_summary <- summary(model) @@ -276,7 +268,6 @@ r2.mlm <- function(model, ...) { } - #' @export r2.glm <- function(model, ci = NULL, verbose = TRUE, ...) { if (!is.null(ci) && !is.na(ci)) { @@ -291,7 +282,7 @@ r2.glm <- function(model, ci = NULL, verbose = TRUE, ...) { if (info$family %in% c("gaussian", "inverse.gaussian")) { out <- r2.default(model, ...) } else if (info$is_logit && info$is_bernoulli) { - out <- list("R2_Tjur" = r2_tjur(model, ...)) + out <- list(R2_Tjur = r2_tjur(model, model_info = info, ...)) attr(out, "model_type") <- "Logistic" names(out$R2_Tjur) <- "Tjur's R2" class(out) <- c("r2_pseudo", class(out)) @@ -301,7 +292,7 @@ r2.glm <- function(model, ci = NULL, verbose = TRUE, ...) { } out <- NULL } else { - out <- list("R2_Nagelkerke" = r2_nagelkerke(model, ...)) + out <- list(R2_Nagelkerke = r2_nagelkerke(model, ...)) names(out$R2_Nagelkerke) <- "Nagelkerke's R2" attr(out, "model_type") <- "Generalized Linear" class(out) <- c("r2_pseudo", class(out)) @@ -313,6 +304,14 @@ r2.glm <- function(model, ci = NULL, verbose = TRUE, ...) { r2.glmx <- r2.glm +#' @export +r2.nestedLogit <- function(model, ci = NULL, verbose = TRUE, ...) { + out <- list(R2_Tjur = r2_tjur(model, ...)) + attr(out, "model_type") <- "Logistic" + class(out) <- c("r2_pseudo", class(out)) + out +} + # mfx models --------------------- @@ -358,7 +357,7 @@ r2.model_fit <- r2.logitmfx #' @export r2.BBreg <- function(model, ...) { - out <- list("R2_CoxSnell" = r2_coxsnell(model)) + out <- list(R2_CoxSnell = r2_coxsnell(model)) names(out$R2_CoxSnell) <- "Cox & Snell's R2" class(out) <- c("r2_pseudo", class(out)) out @@ -378,7 +377,7 @@ r2.bayesx <- r2.BBreg #' @export r2.censReg <- function(model, ...) { - out <- list("R2_Nagelkerke" = r2_nagelkerke(model)) + out <- list(R2_Nagelkerke = r2_nagelkerke(model)) names(out$R2_Nagelkerke) <- "Nagelkerke's R2" class(out) <- c("r2_pseudo", class(out)) out @@ -507,8 +506,8 @@ r2.wbm <- function(model, tolerance = 1e-5, ...) { names(r2_marginal) <- "Marginal R2" out <- list( - "R2_conditional" = r2_conditional, - "R2_marginal" = r2_marginal + R2_conditional = r2_conditional, + R2_marginal = r2_marginal ) attr(out, "model_type") <- "Fixed Effects" @@ -531,8 +530,8 @@ r2.sem <- function(model, ...) { structure( class = "r2_nakagawa", list( - "R2_conditional" = r2_conditional, - "R2_marginal" = r2_marginal + R2_conditional = r2_conditional, + R2_marginal = r2_marginal ) ) } @@ -563,12 +562,10 @@ r2.gam <- function(model, ...) { # gamlss inherits from gam, and summary.gamlss prints results automatically printout <- utils::capture.output(s <- summary(model)) # nolint - if (!is.null(s$r.sq)) { - list( - R2 = c(`Adjusted R2` = s$r.sq) - ) - } else { + if (is.null(s$r.sq)) { NextMethod() + } else { + list(R2 = c(`Adjusted R2` = s$r.sq)) } } @@ -603,7 +600,7 @@ r2.rma <- function(model, ...) { #' @export r2.feis <- function(model, ...) { out <- list( - R2 = c(`R2` = model$r2), + R2 = c(R2 = model$r2), R2_adjusted = c(`adjusted R2` = model$adj.r2) ) @@ -654,7 +651,7 @@ r2.fixest_multi <- function(model, ...) { r2.felm <- function(model, ...) { model_summary <- summary(model) out <- list( - R2 = c(`R2` = model_summary$r2), + R2 = c(R2 = model_summary$r2), R2_adjusted = c(`adjusted R2` = model_summary$r2adj) ) @@ -668,7 +665,7 @@ r2.felm <- function(model, ...) { #' @export r2.iv_robust <- function(model, ...) { out <- list( - R2 = c(`R2` = model$r.squared), + R2 = c(R2 = model$r.squared), R2_adjusted = c(`adjusted R2` = model$adj.r.squared) ) @@ -682,7 +679,7 @@ r2.iv_robust <- function(model, ...) { r2.ivreg <- function(model, ...) { model_summary <- summary(model) out <- list( - R2 = c(`R2` = model_summary$r.squared), + R2 = c(R2 = model_summary$r.squared), R2_adjusted = c(`adjusted R2` = model_summary$adj.r.squared) ) @@ -694,7 +691,7 @@ r2.ivreg <- function(model, ...) { #' @export r2.bigglm <- function(model, ...) { - out <- list("R2_CoxSnell" = summary(model)$rsq) + out <- list(R2_CoxSnell = summary(model)$rsq) names(out$R2_CoxSnell) <- "Cox & Snell's R2" class(out) <- c("r2_pseudo", class(out)) out @@ -728,7 +725,7 @@ r2.biglm <- function(model, ...) { r2.lmrob <- function(model, ...) { model_summary <- summary(model) out <- list( - R2 = c(`R2` = model_summary$r.squared), + R2 = c(R2 = model_summary$r.squared), R2_adjusted = c(`adjusted R2` = model_summary$adj.r.squared) ) @@ -749,10 +746,10 @@ r2.mmclogit <- function(model, ...) { #' @export r2.Arima <- function(model, ...) { - if (!requireNamespace("forecast", quietly = TRUE)) { - list(R2 = NA) - } else { + if (requireNamespace("forecast", quietly = TRUE)) { list(R2 = stats::cor(stats::fitted(model), insight::get_data(model, verbose = FALSE))^2) + } else { + list(R2 = NA) } } @@ -762,8 +759,8 @@ r2.Arima <- function(model, ...) { r2.plm <- function(model, ...) { model_summary <- summary(model) out <- list( - "R2" = c(`R2` = model_summary$r.squared[1]), - "R2_adjusted" = c(`adjusted R2` = model_summary$r.squared[2]) + R2 = c(R2 = model_summary$r.squared[1]), + R2_adjusted = c(`adjusted R2` = model_summary$r.squared[2]) ) attr(out, "model_type") <- "Panel Data" @@ -779,8 +776,8 @@ r2.selection <- function(model, ...) { return(NULL) } out <- list( - "R2" = c(`R2` = model_summary$rSquared$R2), - "R2_adjusted" = c(`adjusted R2` = model_summary$rSquared$R2adj) + R2 = c(R2 = model_summary$rSquared$R2), + R2_adjusted = c(`adjusted R2` = model_summary$rSquared$R2adj) ) attr(out, "model_type") <- "Tobit 2" @@ -795,7 +792,7 @@ r2.svyglm <- function(model, ...) { rsq.adjust <- 1 - ((1 - rsq) * (model$df.null / model$df.residual)) out <- list( - R2 = c(`R2` = rsq), + R2 = c(R2 = rsq), R2_adjusted = c(`adjusted R2` = rsq.adjust) ) @@ -807,7 +804,7 @@ r2.svyglm <- function(model, ...) { #' @export r2.vglm <- function(model, ...) { - out <- list("R2_McKelvey" = r2_mckelvey(model)) + out <- list(R2_McKelvey = r2_mckelvey(model)) names(out$R2_McKelvey) <- "McKelvey's R2" class(out) <- c("r2_pseudo", class(out)) out @@ -820,7 +817,7 @@ r2.vgam <- r2.vglm #' @export r2.DirichletRegModel <- function(model, ...) { - out <- list("R2_Nagelkerke" = r2_nagelkerke(model)) + out <- list(R2_Nagelkerke = r2_nagelkerke(model)) names(out$R2_Nagelkerke) <- "Nagelkerke's R2" class(out) <- c("r2_pseudo", class(out)) out diff --git a/R/r2_coxsnell.R b/R/r2_coxsnell.R index c151136d4..cdb73dea7 100644 --- a/R/r2_coxsnell.R +++ b/R/r2_coxsnell.R @@ -88,6 +88,26 @@ r2_coxsnell.glm <- function(model, verbose = TRUE, ...) { #' @export r2_coxsnell.BBreg <- r2_coxsnell.glm + +#' @export +r2_coxsnell.nestedLogit <- function(model, ...) { + n <- insight::n_obs(model, disaggregate = TRUE) + stats::setNames( + lapply(names(model$models), function(i) { + m <- model$models[[i]] + # if no deviance, return NA + if (is.null(m$deviance)) { + return(NA) + } + r2_coxsnell <- (1 - exp((m$deviance - m$null.deviance) / n[[i]])) + names(r2_coxsnell) <- "Cox & Snell's R2" + r2_coxsnell + }), + names(model$models) + ) +} + + #' @export r2_coxsnell.mclogit <- function(model, ...) { insight::check_if_installed("mclogit", reason = "to calculate R2") diff --git a/R/r2_nagelkerke.R b/R/r2_nagelkerke.R index bb6230f22..85bcc6e8f 100644 --- a/R/r2_nagelkerke.R +++ b/R/r2_nagelkerke.R @@ -77,6 +77,26 @@ r2_nagelkerke.glm <- function(model, verbose = TRUE, ...) { #' @export r2_nagelkerke.BBreg <- r2_nagelkerke.glm + +#' @export +r2_nagelkerke.nestedLogit <- function(model, ...) { + n <- insight::n_obs(model, disaggregate = TRUE) + stats::setNames( + lapply(names(model$models), function(i) { + m <- model$models[[i]] + # if no deviance, return NA + if (is.null(m$deviance)) { + return(NA) + } + r2_nagelkerke <- (1 - exp((m$deviance - m$null.deviance) / n[[i]])) / (1 - exp(-m$null.deviance / n[[i]])) + names(r2_nagelkerke) <- "Nagelkerke's R2" + r2_nagelkerke + }), + names(model$models) + ) +} + + #' @export r2_nagelkerke.bife <- function(model, ...) { r2_nagelkerke <- r2_coxsnell(model) / (1 - exp(-model$null_deviance / insight::n_obs(model))) diff --git a/R/r2_tjur.R b/R/r2_tjur.R index e0aa6117c..6a044e299 100644 --- a/R/r2_tjur.R +++ b/R/r2_tjur.R @@ -23,6 +23,11 @@ #' #' @export r2_tjur <- function(model, ...) { + UseMethod("r2_tjur") +} + +#' @export +r2_tjur.default <- function(model, ...) { info <- list(...)$model_info if (is.null(info)) { info <- suppressWarnings(insight::model_info(model, verbose = FALSE)) @@ -50,3 +55,29 @@ r2_tjur <- function(model, ...) { names(tjur_d) <- "Tjur's R2" tjur_d } + +#' @export +r2_tjur.nestedLogit <- function(model, ...) { + resp <- insight::get_response(model, dichotomies = TRUE, verbose = FALSE) + + stats::setNames( + lapply(names(model$models), function(i) { + y <- resp[[i]] + m <- model$models[[i]] + pred <- stats::predict(m, type = "response") + # delete pred for cases with missing residuals + if (anyNA(stats::residuals(m))) { + pred <- pred[!is.na(stats::residuals(m))] + } + categories <- unique(y) + m1 <- mean(pred[which(y == categories[1])], na.rm = TRUE) + m2 <- mean(pred[which(y == categories[2])], na.rm = TRUE) + + tjur_d <- abs(m2 - m1) + + names(tjur_d) <- "Tjur's R2" + tjur_d + }), + names(model$models) + ) +} diff --git a/man/check_itemscale.Rd b/man/check_itemscale.Rd index 7fa487ab5..7f790b1d2 100644 --- a/man/check_itemscale.Rd +++ b/man/check_itemscale.Rd @@ -61,7 +61,5 @@ check_itemscale(pca) \item Briggs SR, Cheek JM (1986) The role of factor analysis in the development and evaluation of personality scales. Journal of Personality, 54(1), 106-148. doi: 10.1111/j.1467-6494.1986.tb00391.x -\item Trochim WMK (2008) Types of Reliability. -(\href{https://conjointly.com/kb/types-of-reliability/}{web}) } } diff --git a/tests/testthat/_snaps/check_distribution.md b/tests/testthat/_snaps/check_distribution.md new file mode 100644 index 000000000..11187a23b --- /dev/null +++ b/tests/testthat/_snaps/check_distribution.md @@ -0,0 +1,21 @@ +# check_distribution + + Code + print(out) + Output + # Distribution of Model Family + + Predicted Distribution of Residuals + + Distribution Probability + cauchy 94% + lognormal 3% + weibull 3% + + Predicted Distribution of Response + + Distribution Probability + lognormal 47% + gamma 44% + beta-binomial 3% + diff --git a/tests/testthat/_snaps/nestedLogit.md b/tests/testthat/_snaps/nestedLogit.md new file mode 100644 index 000000000..f5c9a0bdd --- /dev/null +++ b/tests/testthat/_snaps/nestedLogit.md @@ -0,0 +1,12 @@ +# model_performance + + Code + model_performance(mnl) + Output + # Indices of model performance + + Response | AIC | BIC | RMSE | Sigma | R2 + ---------------------------------------------------- + work | 325.733 | 336.449 | 0.456 | 1.000 | 0.138 + full | 110.495 | 118.541 | 0.398 | 1.000 | 0.333 + diff --git a/tests/testthat/test-binned_residuals.R b/tests/testthat/test-binned_residuals.R new file mode 100644 index 000000000..4aa69e0ec --- /dev/null +++ b/tests/testthat/test-binned_residuals.R @@ -0,0 +1,67 @@ +test_that("binned_residuals", { + data(mtcars) + model <- glm(vs ~ wt + mpg, data = mtcars, family = "binomial") + result <- binned_residuals(model) + expect_named( + result, + c("xbar", "ybar", "n", "x.lo", "x.hi", "se", "ci_range", "CI_low", "CI_high", "group") + ) + expect_equal( + result$xbar, + c(0.03786, 0.09514, 0.25911, 0.47955, 0.71109, 0.97119), + tolerance = 1e-4 + ) + expect_equal( + result$ybar, + c(-0.03786, -0.09514, 0.07423, -0.07955, 0.28891, -0.13786), + tolerance = 1e-4 + ) +}) + + +test_that("binned_residuals, n_bins", { + data(mtcars) + model <- glm(vs ~ wt + mpg, data = mtcars, family = "binomial") + result <- binned_residuals(model, n_bins = 10) + expect_named( + result, + c("xbar", "ybar", "n", "x.lo", "x.hi", "se", "ci_range", "CI_low", "CI_high", "group") + ) + expect_equal( + result$xbar, + c( + 0.02373, 0.06301, 0.08441, 0.17907, 0.29225, 0.44073, 0.54951, + 0.69701, 0.9168, 0.99204 + ), + tolerance = 1e-4 + ) + expect_equal( + result$ybar, + c( + -0.02373, -0.06301, -0.08441, -0.17907, 0.20775, -0.1074, 0.11715, + 0.30299, -0.25014, 0.00796 + ), + tolerance = 1e-4 + ) +}) + + +test_that("binned_residuals, terms", { + data(mtcars) + model <- glm(vs ~ wt + mpg, data = mtcars, family = "binomial") + result <- binned_residuals(model, term = "mpg") + expect_named( + result, + c("xbar", "ybar", "n", "x.lo", "x.hi", "se", "ci_range", "CI_low", "CI_high", "group") + ) + expect_equal( + result$xbar, + c(12.62, 15.34, 18.1, 20.9, 22.875, 30.06667), + tolerance = 1e-4 + ) + expect_equal( + result$ybar, + c(-0.05435, -0.07866, 0.13925, -0.11861, 0.27763, -0.13786), + tolerance = 1e-4 + ) +}) diff --git a/tests/testthat/test-check_autocorrelation.R b/tests/testthat/test-check_autocorrelation.R new file mode 100644 index 000000000..b389f3985 --- /dev/null +++ b/tests/testthat/test-check_autocorrelation.R @@ -0,0 +1,7 @@ +test_that("check_autocorrelation", { + data(mtcars) + m <- lm(mpg ~ wt + cyl + gear + disp, data = mtcars) + set.seed(123) + out <- check_autocorrelation(m) + expect_equal(as.vector(out), 0.316, ignor_attr = TRUE, tolerance = 1e-2) +}) diff --git a/tests/testthat/test-check_distribution.R b/tests/testthat/test-check_distribution.R new file mode 100644 index 000000000..e8ab2835f --- /dev/null +++ b/tests/testthat/test-check_distribution.R @@ -0,0 +1,36 @@ +test_that("check_distribution", { + skip_if_not_installed("lme4") + skip_if_not_installed("randomForest") + data(sleepstudy, package = "lme4") + model <<- lme4::lmer(Reaction ~ Days + (Days | Subject), sleepstudy) + out <- check_distribution(model) + + expect_identical( + out$Distribution, + c( + "bernoulli", "beta", "beta-binomial", "binomial", "cauchy", + "chi", "exponential", "F", "gamma", "half-cauchy", "inverse-gamma", + "lognormal", "neg. binomial (zero-infl.)", "negative binomial", + "normal", "pareto", "poisson", "poisson (zero-infl.)", "tweedie", + "uniform", "weibull" + ) + ) + expect_equal( + out$p_Residuals, + c( + 0, 0, 0, 0, 0.9375, 0, 0, 0, 0, 0, 0, 0.03125, 0, 0, 0, 0, + 0, 0, 0, 0, 0.03125 + ), + tolerance = 1e-4 + ) + expect_equal( + out$p_Response, + c( + 0, 0, 0.03125, 0, 0, 0, 0, 0, 0.4375, 0.03125, 0, 0.46875, + 0.03125, 0, 0, 0, 0, 0, 0, 0, 0 + ), + tolerance = 1e-4 + ) + + expect_snapshot(print(out)) +}) diff --git a/tests/testthat/test-check_heterogeneity_bias.R b/tests/testthat/test-check_heterogeneity_bias.R new file mode 100644 index 000000000..7abc6af30 --- /dev/null +++ b/tests/testthat/test-check_heterogeneity_bias.R @@ -0,0 +1,33 @@ +test_that("check_heterogeneity_bias", { + data(iris) + set.seed(123) + iris$ID <- sample(1:4, nrow(iris), replace = TRUE) # fake-ID + out <- check_heterogeneity_bias(iris, select = c("Sepal.Length", "Petal.Length"), group = "ID") + expect_equal(out, c("Sepal.Length", "Petal.Length"), ignore_attr = TRUE) + expect_output(print(out), "Possible heterogeneity bias due to following predictors: Sepal\\.Length, Petal\\.Length") + + out <- check_heterogeneity_bias(iris, select = ~ Sepal.Length + Petal.Length, group = ~ID) + expect_equal(out, c("Sepal.Length", "Petal.Length"), ignore_attr = TRUE) + expect_output(print(out), "Possible heterogeneity bias due to following predictors: Sepal\\.Length, Petal\\.Length") + + m <- lm(Sepal.Length ~ Petal.Length + Petal.Width + Species + ID, data = iris) + expect_error( + check_heterogeneity_bias(m, select = c("Sepal.Length", "Petal.Length"), group = "ID"), + regex = "no mixed model" + ) + + skip_if_not_installed("lme4") + m <- lme4::lmer(Sepal.Length ~ Petal.Length + Petal.Width + Species + (1 | ID), data = iris) + out <- check_heterogeneity_bias(m, select = c("Sepal.Length", "Petal.Length"), group = "ID") + expect_equal(out, c("Petal.Length", "Petal.Width", "Species"), ignore_attr = TRUE) + expect_output( + print(out), + "Possible heterogeneity bias due to following predictors: Petal\\.Length, Petal\\.Width, Species" + ) + out <- check_heterogeneity_bias(m, select = ~ Sepal.Length + Petal.Length, group = ~ID) + expect_equal(out, c("Petal.Length", "Petal.Width", "Species"), ignore_attr = TRUE) + expect_output( + print(out), + "Possible heterogeneity bias due to following predictors: Petal\\.Length, Petal\\.Width, Species" + ) +}) diff --git a/tests/testthat/test-nestedLogit.R b/tests/testthat/test-nestedLogit.R new file mode 100644 index 000000000..47a852b04 --- /dev/null +++ b/tests/testthat/test-nestedLogit.R @@ -0,0 +1,66 @@ +skip_on_os(c("mac", "linux")) +skip_if(packageVersion("insight") <= "0.19.5.10") +skip_if_not_installed("nestedLogit") +skip_if_not_installed("carData") + +data("Womenlf", package = "carData") +comparisons <- nestedLogit::logits( + work = nestedLogit::dichotomy("not.work", working = c("parttime", "fulltime")), + full = nestedLogit::dichotomy("parttime", "fulltime") +) +mnl <- nestedLogit::nestedLogit( + partic ~ hincome + children, + dichotomies = comparisons, + data = Womenlf +) + +test_that("r2", { + out <- r2(mnl) + expect_equal( + out, + list(R2_Tjur = list( + work = c(`Tjur's R2` = 0.137759452521642), + full = c(`Tjur's R2` = 0.332536937208286) + )), + ignore_attr = TRUE, + tolerance = 1e-4 + ) + + out <- r2_tjur(mnl) + expect_equal( + out, + list( + work = c(`Tjur's R2` = 0.137759452521642), + full = c(`Tjur's R2` = 0.332536937208286) + ), + ignore_attr = TRUE, + tolerance = 1e-4 + ) + + out <- r2_coxsnell(mnl) + expect_equal( + out, + list( + work = c(`Cox & Snell's R2` = 0.129313084315599), + full = c(`Cox & Snell's R2` = 0.308541455410686) + ), + ignore_attr = TRUE, + tolerance = 1e-4 + ) + + out <- r2_nagelkerke(mnl) + expect_equal( + out, + list( + work = c(`Nagelkerke's R2` = 0.174313365512442), + full = c(`Nagelkerke's R2` = 0.418511411473948) + ), + ignore_attr = TRUE, + tolerance = 1e-4 + ) +}) + + +test_that("model_performance", { + expect_snapshot(model_performance(mnl)) +})