R/get.R

In brokenstick: Broken Stick Model for Irregular Longitudinal Data

#' Obtain the knots from a broken stick model
#'
#' @param object An object of class \code{brokenstick}
#' @param whatknots Deprecated. Use `hide` instead.
#' @param what Deprecated. Use `hide` instead.
#' @inheritParams brokenstick
#' @return A vector with knot locations, either both, internal only or
#' boundary only, depending on `hide`.
#' The result is \code{NULL} if \code{object} does not
#' have proper class. Returns \code{numeric(0)} if
#' there are no internal knots.
#' @examples
#' get_knots(fit_200, hide = "bo")
#' @export
get_knots <- function(
  object,
  hide = c("right", "left", "boundary", "internal", "none"),
  whatknots = "all",
  what = "all"
) {
  stopifnot(inherits(object, c("brokenstick")))

  if (!missing(what) || !missing(whatknots)) {
    warning(
      "arguments 'what' and 'whatknots' in 'get_knots()' are deprecated; please use 'hide' instead.",
      call. = FALSE
    )
    whatknots <- what
    object$hide <- switch(
      whatknots,
      all = "none",
      internal = "boundary",
      boundary = "internal",
      dropfirst = "left",
      droplast = "right",
      "right"
    )
  }

  if (!missing(hide)) {
    hide <- match.arg(hide)
  } else {
    hide <- ifelse(is.null(object$hide), "right", object$hide)
  }

  internal <- object$internal
  # legacy for objects created before v2.0
  if (is.null(internal)) {
    internal <- object$knots
  }
  boundary <- object$boundary
  internal <- internal[internal > boundary[1L] & internal < boundary[2L]]

  result <- switch(
    hide,
    none = c(boundary[1L], internal, boundary[2L]),
    internal = boundary,
    boundary = internal,
    left = c(internal, boundary[2L]),
    right = c(boundary[1L], internal)
  )
  return(result)
}

#' Extract Variance and Correlation Components
#'
#' Extracts variance-covariance or correlation matrix from a
#' `brokenstick` object.
#'
#' @param x      Object of class `brokenstick`
#' @param cor    Logical. Should the function return the correlation matrix
#' instead of the covariance matrix? The default is `FALSE`.
#' @param what  Deprecated.
#' @param whatknots Deprecated.
#' @inheritParams brokenstick
#' @return A numeric matrix, possibly with zero rows and columns if no names match
#' @examples
#' f1 <- brokenstick(hgt_z ~ age | id, smocc_200[1:1000, ], knots = 0:2, seed = 1)
#' get_omega(f1, cor = TRUE, hide = "boundary")
#' @export
get_omega <- function(
  x,
  hide = c("right", "left", "boundary", "internal", "none"),
  cor = FALSE,
  whatknots = "all",
  what = "cov"
) {
  stopifnot(inherits(x, "brokenstick"))

  if (!missing(what)) {
    warning(
      "argument 'what' in 'get_omega()' is deprecated; please use 'cor' instead.",
      call. = FALSE
    )
    if (missing(cor)) {
      cor <- ifelse(what == "cor", TRUE, FALSE)
    }
  }

  if (!missing(whatknots)) {
    warning(
      "argument 'whatknots' in 'get_omega()' is deprecated; please use 'hide' instead.",
      call. = FALSE
    )
    whatknots <- what
    x$hide <- switch(
      whatknots,
      all = "none",
      internal = "boundary",
      boundary = "internal",
      dropfirst = "left",
      droplast = "right",
      "right"
    )
  }

  if (!missing(hide)) {
    hide <- match.arg(hide)
  } else {
    hide <- ifelse(is.null(x$hide), "right", x$hide)
  }

  omega <- x$omega
  v <- colnames(omega)
  nameset <- switch(
    hide,
    none = v,
    boundary = v[c(-1L, -length(v))],
    internal = v[c(1L, length(v))],
    left = v[-1L],
    right = v[-length(v)]
  )
  if (length(nameset)) {
    omega <- omega[nameset, nameset, drop = FALSE]
  } else {
    omega <- matrix(NA_real_, 0L, 0L)
  }

  # return correlation matrix if requested
  if (cor) {
    if (dim(omega)[1L]) {
      return(cov2cor(omega))
    }
  }

  # covariance matrix
  return(omega)
}


#' Obtain proportion of explained variance from a broken stick model
#'
#' @param object An object of class \code{brokenstick}
#' @param newdata Data on which `r.squared` must be calculated
#' @return Proportion of explained variance
#' @examples
#' get_r2(fit_200)
#' get_r2(fit_200_light, newdata = smocc_200)
#' @export
get_r2 <- function(object, newdata = NULL) {
  newdata <- get_newdata(object, newdata)
  p <- predict(object, newdata = newdata)
  nd <- newdata %>%
    select(object$names$y) %>%
    bind_cols(p) %>%
    tidyr::drop_na()
  return(cor(nd[[".pred"]], nd[[object$names$y]])^2)
}


get_newdata <- function(x, newdata) {
  # sets the newdata argument
  stopifnot(inherits(x, "brokenstick"))
  if (is.null(newdata) && x$light) {
    stop(
      "Argument 'newdata' is required for a light brokenstick object.",
      call. = FALSE
    )
  }
  if (is.null(newdata) && !x$light) {
    newdata <- x$data
  }
  stopifnot(is.data.frame(newdata) || is.matrix(newdata))
  return(newdata)
}