File: select_parameters.R

package info (click to toggle)
r-cran-parameters 0.24.2-2
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 3,852 kB
  • sloc: sh: 16; makefile: 2
file content (119 lines) | stat: -rw-r--r-- 3,681 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
#' Automated selection of model parameters
#'
#' This function performs an automated selection of the 'best' parameters,
#' updating and returning the "best" model.
#'
#' @param model A statistical model (of class `lm`, `glm`, or `merMod`).
#' @param ... Arguments passed to or from other methods.
#'
#' @section Classical lm and glm:
#' For frequentist GLMs, `select_parameters()` performs an AIC-based stepwise
#' selection.
#'
#' @section Mixed models:
#' For mixed-effects models of class `merMod`, stepwise selection is based on
#' [`cAIC4::stepcAIC()`]. This step function only searches the "best" model
#' based on the random-effects structure, i.e. `select_parameters()` adds or
#' excludes random-effects until the cAIC can't be improved further.
#'
#' @examplesIf requireNamespace("lme4")
#' model <- lm(mpg ~ ., data = mtcars)
#' select_parameters(model)
#'
#' model <- lm(mpg ~ cyl * disp * hp * wt, data = mtcars)
#' select_parameters(model)
#' \donttest{
#' # lme4 -------------------------------------------
#' model <- lme4::lmer(
#'   Sepal.Width ~ Sepal.Length * Petal.Width * Petal.Length + (1 | Species),
#'   data = iris
#' )
#' select_parameters(model)
#' }
#'
#' @return The model refitted with optimal number of parameters.
#' @export
select_parameters <- function(model, ...) {
  UseMethod("select_parameters")
}


#' @rdname select_parameters
#' @param k The multiple of the number of degrees of freedom used for the penalty.
#' Only `k = 2` gives the genuine AIC: `k = log(n)` is sometimes referred to as
#' BIC or SBC.
#' @inheritParams stats::step
#' @export
select_parameters.lm <- function(model,
                                 direction = "both",
                                 steps = 1000,
                                 k = 2,
                                 ...) {
  junk <- utils::capture.output({
    best <- stats::step(
      model,
      trace = 0,
      direction = direction,
      steps = steps,
      k = k,
      ...
    )
  })

  best
}


#' @rdname select_parameters
#' @export
select_parameters.merMod <- function(model,
                                     direction = "backward",
                                     steps = 1000,
                                     ...) {
  insight::check_if_installed("cAIC4")

  # Find slope and group candidates
  # data <- insight::get_data(model)
  # factors <- names(data[sapply(data, is.factor)])
  # if(length(factors) == 0){
  #   factors <- NULL
  # }
  # nums <- names(data[sapply(data, is.numeric)])
  # if(length(nums) == 0){
  #   nums <- NULL
  # }

  factors <- unique(c(
    insight::find_random(model, split_nested = FALSE, flatten = TRUE),
    insight::find_random(model, split_nested = TRUE, flatten = TRUE)
  ))
  factors <- gsub(":", "/", factors, fixed = TRUE)


  best <- suppressMessages(
    suppressWarnings(
      cAIC4::stepcAIC(
        model,
        # slopeCandidates = nums,
        groupCandidates = factors,
        direction = direction,
        steps = steps,
        allowUseAcross = TRUE
      )$finalModel
    )
  )


  # Using MuMIn's dredge(): works nicely BUT throws unnecessary warnings and
  # requires to set global options for na.action even tho no NaNs.
  # The code is here: https://github.com/cran/MuMIn/blob/master/R/dredge.R Maybe it could be reimplemented?
  # insight::check_if_installed("MuMIn")
  # model <- lmer(
  #  Sepal.Width ~ Sepal.Length * Petal.Width * Petal.Length + (1 | Species),
  #  data = iris,
  #  na.action = na.fail
  # )
  # summary(MuMIn::get.models(MuMIn::dredge(model), 1)[[1]])

  best
}