R/PLR-predict.R

In LorenzRegression: Lorenz and Penalized Lorenz Regressions

#' Prediction and fitted values for the penalized Lorenz regression
#'
#' \code{predict} provides predictions for an object of class \code{"PLR"},
#' while \code{fitted} extracts the fitted values.
#'
#' @aliases predict.PLR_boot predict.PLR_cv fitted.PLR fitted.PLR_boot fitted.PLR_cv
#' @param object An object of S3 class \code{"PLR"}. The object might also have S3 classes \code{"PLR_boot"} and/or \code{"PLR_cv"} (both inherit from class \code{"PLR"})
#' @param newdata An optional data frame in which to look for variables with which to predict. If omitted, the original data are used.
#' @param type A character string indicating the type of prediction or fitted values. Possible values are \code{"index"} (the default) and response.
#' In the first case, only the index of the single-index model is estimated.
#' In the second case, the "full" conditional expectation of the response given the covariates is estimated.
#' @param pars.idx What grid and penalty parameters should be used for parameter selection. Either a character string specifying the selection method, where the possible values are:
#' \itemize{
#'    \item \code{"BIC"} (default) - Always available.
#'    \item \code{"Boot"} - Available if \code{object} inherits from \code{"PLR_boot"}.
#'    \item \code{"CV"} - Available if \code{object} inherits from \code{"PLR_cv"}.
#' }
#' Or a numeric vector of length 2, where the first element is the index of the grid parameter and the second is the index of the penalty parameter.
#' @param ... Additional arguments passed to the function \code{\link{Rearrangement.estimation}}.
#'
#' @return A vector of predictions for \code{predict}, or a vector of fitted values for \code{fitted}.
#'
#' @details
#' The \code{type} argument distinguishes between two types of prediction outputs, aligned with the goals of the penalized Lorenz regression.
#' When \code{type = "index"}, the function returns the estimated index \eqn{X^\top \theta} of the single-index model. This index captures the full ordering structure of the conditional expectation and is sufficient for computing the explained Gini coefficient, which is the primary focus of the method. Crucially, this estimation does not require recovering the full nonparametric link function.
#' When \code{type = "response"}, the function estimates the full conditional expectation \eqn{\mathbb{E}[Y | X]} by performing a second-stage estimation of the link function via \code{\link{Rearrangement.estimation}}. This is useful if fitted or predicted response values are needed for other purposes.
#'
#' @seealso \code{\link{Lorenz.Reg}}, \code{\link{Rearrangement.estimation}}
#'
#' @examples
#' ## For examples see example(Lorenz.Reg), example(Lorenz.boot) and example(PLR.CV)
#'
#' @importFrom stats terms delete.response model.frame
#'
#' @method predict PLR
#' @export

predict.PLR <- function(object, newdata, type=c("index","response"), pars.idx = "BIC", ...){

  type <- match.arg(type)

  if((is.numeric(pars.idx) & length(pars.idx)==2)){
    lth1 <- length(object$path)
    if(pars.idx[1] > lth1) stop("Index of grid parameter is out of bond.")
    lth2 <- ncol(object$path[[pars.idx[1]]])
    if(pars.idx[2] > lth2) stop("Index of lambda parameter is out of bond.")
  }else if(pars.idx == "BIC"){
    pars.idx <- c(object$grid.idx["BIC"],object$lambda.idx["BIC"])
  }else if(pars.idx == "Boot"){
    stop("object is not of class 'PLR_boot'. Therefore pars.idx cannot be set to 'Boot'.")
  }else if(pars.idx == "CV"){
    stop("object is not of class 'PLR_cv'. Therefore pars.idx cannot be set to 'CV'.")
  }else{
    stop("pars.idx does not have the correct format")
  }

  predict_PLR(object, newdata, type, pars.idx, ...)

}

#' @method predict PLR_boot
#' @export

predict.PLR_boot <- function(object, newdata, type=c("index","response"), pars.idx = "BIC", ...){

  type <- match.arg(type)

  if(all(pars.idx == "Boot")) pars.idx <- c(object$grid.idx["Boot"],object$lambda.idx["Boot"])
  NextMethod("predict")

}

#' @method predict PLR_cv
#' @export

predict.PLR_cv <- function(object, newdata, type=c("index","response"), pars.idx = "BIC", ...){

  type <- match.arg(type)

  if(all(pars.idx == "CV")) pars.idx <- c(object$grid.idx["CV"],object$lambda.idx["CV"])
  NextMethod("predict")

}

predict_PLR <- function(object, newdata, type, pars.idx, ...){

  # Data (re)-construction
  tt <- terms(object)
  noData <- (missing(newdata) || is.null(newdata))
  if(noData){
    x <- object$x
  }else{
    Terms <- delete.response(tt)
    m <- model.frame(Terms, newdata, xlev = object$xlevels)
    x <- model_matrix_PLR(Terms,m)
  }

  # Retrieving theta and index
  l <- ncol(object$x)
  pth <- object$path[[pars.idx[1]]][,pars.idx[2]]
  object$theta <- pth[(length(pth)-l+1):length(pth)]
  object$index <- as.vector(object$theta%*%t(object$x)) # Necessarily on the original x
  index <- object$theta%*%t(x) # on the x used for predictions

  # Defining the predictor
  if(type=="index"){
    predictor <- as.vector(index)
  }else{
    predictor <- Rearrangement.estimation(object$y, object$index, t=as.vector(index), weights=object$weights, ...)$H
    names(predictor) <- NULL
  }

  return(predictor)

}

#' @method fitted PLR
#' @rdname predict.PLR
#' @export

fitted.PLR <- function(object, type=c("index","response"), pars.idx = "BIC", ...){

  predict.PLR(object, type = type, pars.idx = pars.idx, ...)

}

#' @method fitted PLR_boot
#' @export

fitted.PLR_boot <- function(object, type=c("index","response"), pars.idx = "BIC", ...){

  if(all(pars.idx == "Boot")) pars.idx <- c(object$grid.idx["Boot"],object$lambda.idx["Boot"])
  NextMethod("fitted")

}

#' @method fitted PLR_cv
#' @export

fitted.PLR_cv <- function(object, type=c("index","response"), pars.idx = "BIC", ...){

  if(all(pars.idx == "CV")) pars.idx <- c(object$grid.idx["CV"],object$lambda.idx["CV"])
  NextMethod("fitted")

}

#' Residuals for the penalized Lorenz regression
#'
#' \code{residuals} provides residuals for an object of class \code{"PLR"}.
#'
#' @aliases residuals.PLR_boot residuals.PLR_cv
#' @param object An object of class \code{"PLR"}.
#' @param pars.idx What grid and penalty parameters should be used for parameter selection. Either a character string specifying the selection method, where the possible values are:
#' \itemize{
#'    \item \code{"BIC"} (default) - Always available.
#'    \item \code{"Boot"} - Available if \code{object} inherits from \code{"PLR_boot"}.
#'    \item \code{"CV"} - Available if \code{object} inherits from \code{"PLR_cv"}.
#' }
#' Or a numeric vector of length 2, where the first element is the index of the grid parameter and the second is the index of the penalty parameter.
#' @param ... Additional arguments passed to the function \code{\link{Rearrangement.estimation}}.
#'
#' @return A vector of residuals.
#'
#' @details Computing residuals entail to estimate the link function of the single-index model. This is done via the function \code{\link{Rearrangement.estimation}}.
#'
#' @seealso \code{\link{Lorenz.Reg}}, \code{\link{Rearrangement.estimation}}
#'
#' @examples
#' ## For examples see example(Lorenz.Reg), example(Lorenz.boot) and example(PLR.CV)
#'
#' @method residuals PLR
#' @export

residuals.PLR <- function(object, pars.idx = "BIC", ...){

  yhat <- fitted.PLR(object, type = "response", pars.idx = pars.idx, ...)
  y <- object$y
  r <- y - yhat
  return(r)

}

#' @method residuals PLR_boot
#' @export

residuals.PLR_boot <- function(object, pars.idx = "BIC", ...){
  if(all(pars.idx == "Boot")) pars.idx <- c(object$grid.idx["Boot"],object$lambda.idx["Boot"])
  NextMethod("residuals")
}

#' @method residuals PLR_cv
#' @export

residuals.PLR_cv <- function(object, pars.idx = "BIC", ...){
  if(all(pars.idx == "CV")) pars.idx <- c(object$grid.idx["CV"],object$lambda.idx["CV"])
  NextMethod("residuals")
}