R/portfolio_selection.R

In benSepanski/PortfolioSelection: Test and develop portfolio selection algorithms

# PortfolioSelector -------------------------------------------------------


#' An S4 class to represent a portfolio selection strategy.
#'
#' @slot price_relatives A data frame where each row represents a
#'     trading period, and each column represents a stock.
#'     Each entry is a price relative: if the price
#'     of a stock is \eqn{p_t} at time period
#'     \eqn{t}, it has price relative \eqn{\frac{p_t}{p_{t-1}}}
#'     at time priod \eqn{t}.
#' @slot portfolio A vector of entries \eqn{b_i} where
#'     \eqn{b_i} represents the proportion of your portfolio's
#'     total wealth which is allocated to asset \eqn{i}
#'     in the current trading period
#' @slot trading_period An row index into \code{@price_relatives}
#'     representing the current trading period.
#'
methods::setClass("PortfolioSelector",
        slots = c(price_relatives = "data.frame",
                  portfolio = "numeric",
                  trading_period = "numeric",
                  transaction_rate = "numeric")
)

#' @describeIn Verifies a PortfolioSelector object
#'
#' @title check_PortfolioSelector
#'
#' Makes sure that a PortfolioSelector object has an
#' existent price_relatives data frame, a trading_period
#' which is a valid row in the data frame, and
#' a portfolio whose length matches the number of columns
#' of the data frame. Also makes sure that
#' the portfolio is non-negative with entries that sum to 1j
#'
#' @param object the PortfolioSelector object to verify
#' @return \code{TRUE} if no errors, otherwise a vector of strings
#'     describing the errors.
methods::setValidity("PortfolioSelector", function(object) {
  errors <- character()
  # make sure transaction rate is a single numeric in [0,1]
  tr <- object@transaction_rate
  if(any(is.na(tr))) {
    errors <- c(errors,
                glue::glue("transaction_rate has missing values."))
  }
  else if(!rlang::is_scalar_double(tr)) {
    errors <- c(errors,
                glue::glue("transaction_rate must be a scalar double."))
  }
  else if(tr < 0 || tr > 1) {
    errors <- c(errors,
                glue::glue("transaction_rate of {tr} is not in the unit \\
interval [0,1]."))
  }
  # Make sure the price relatives are not missing
  ntime_periods <- nrow(object@price_relatives)
  nassets <- ncol(object@price_relatives)
  if(ntime_periods <= 0 || nassets <= 0) {
    msg <- glue::glue("price_relatives has {ntime_periods} rows \\
and {nassets} columns, but must have at least 1 row and at least \\
1 column.")
    errors <- c(errors, msg)
  }
  # make sure the trading period is a single whole number which
  # is not missing and is valid
  tp <- object@trading_period
  if(any(is.na(tp))) {
    errors <- c(errors, glue::glue("trading_period has missing values."))
  }
  else if(length(tp) != 1) {
    errors <- c(errors, glue::glue("trading_period must be a scalar value."))
  }
  else if(!is_whole_number(tp)) {
    errors <- c(errors, glue::glue("trading_period must be (close to) a \\
whole number."))
  }
  else if(tp > ntime_periods || tp < 0) {
    errors <- c(errors, glue::glue("trading_period out of bounds: must be in \\
range [0, {ntime_periods}]."))
  }
  # Make sure the portfolio is not missing and is the right shape
  if(any(is.na(object@portfolio))) {
    msg <- "portfolio has missing entries."
    errors <- c(errors, msg)
  }
  if(length(object@portfolio) != nassets) {
    msg <- glue::glue("portfolio is length {length(object@portfolio)}, which \\
does not equal the number of columns in price_relatives ({nassets}).")
    errors <- c(errors, msg)
  }
  # make sure portfolio is a partition of 1
  if(any(object@portfolio < 0)) {
    msg <- "portfolio has negative entries."
    errors <- c(errors, msg)
  }
  if(abs(sum(object@portfolio) - 1) > 1e-7) {
    msg <- "portfolio sums to {sum(object@portfolio)}, but should sum to 1."
    errors <- c(errors, msg)
  }

  if(length(errors) == 0) TRUE else errors
})


#' Return the next portfolio according to the selector's strategy.
#'
#' During \code{pf_select}'s current trading period, trade
#' according to \code{pf_select}'s strategy. Namely,
#' using the prices at trading period \code{pf_select@@trading_period},
#' return \code{pf_select}'s recommended new portfolio.
#'
#' @param pf_select A PortfolioSelector object
#'     used to select the next portfolio
#' @return Returns the next portfolio according to the selector,
#'     with the portfolio distributed according to the prices
#'     at trading period \code{pf_select@@trading_period}.
#'
methods::setGeneric("recommend_next_portfolio", function(pf_selector) {
  standardGeneric("recommend_next_portfolio")
})


# BuyAndHold --------------------------------------------------------------


#' A \code{PortfolioSelector} class which never trades stocks
#'
#' This strategy takes its initial portfolio and does not
#' sell or buy stocks. The portfolio only changes according
#' to price adjustments.
#'
methods::setClass("BuyAndHold", contains = "PortfolioSelector")

#' @describeIn BuyAndHold Adjusts portfolio to next trading period prices
#'
#' Returns portfolio for the next trading period by not
#' selling or buying any stocks, simply adjusting the prices
#'
methods::setMethod("recommend_next_portfolio", signature = "BuyAndHold",
function(pf_selector) {
  pf_selector@portfolio
})


# LOAD --------------------------------------------------------------------


#' Implements the LOAD Online PS System
#'
#' A class to run the LOAD online PS System: For each stock
#' regresses the prices in the last time_window many trading
#' periods and thresholds the slope to classify the stock as
#' with momentum or without momentum. It uses this classification
#' to predict the price vector, then chooses the closest portfolio
#' which is predicted to return at least wealth_factor_threshold
#' over the next period
#'
#' This method was developed in the following paper:
#' @seealso \url{https://www.sciencedirect.com/science/article/abs/pii/S0950705119303922#b19}
#'
#' @slot decay_factor \eqn{\alpha} in the referenced paper, LOAD predicts
#'     that stocks regressing to the mean have approximate price
#'     \eqn{MA_t = \alpha p_t + (1-\alpha)MA_{t-1}}, \eqn{MA_1 = p_1}
#' @slot regularization_factor \eqn{\lambda} in the referenced paper,
#'     the regularization coeffecient for weight decay when
#'     regressing the prices in the recent time window.
#'     Must be non-negative
#' @slot time_window \eqn{w} in the referenced paper.
#'     The number of trading periods (including
#'     the most recent) to regress on (must be at least 2)
#' @slot momentum_threshold \eqn{\eta} in the referenced paper.
#'     If the regressed slope is greater than \eqn{\eta} we say the
#'     stock has momentum. This must be greater than 0.
#' @slot wealth_factor_threshold \eqn{\epsilon} in the referenced
#'     paper. Once a price relative for the next
#'     trading period is predicted, this is the minimum return
#'     our traded portfolio should satisfy. This must be greater
#'     than 0.
#' @slot prices (OPTIONAL--THROWS WARNING IF NOT SUPPLIED)
#'     The prices during the trading period. If not supplied,
#'     it is assumed all the initial prices are 1.
#' @slot price_means (OPTIONAL)
#'     The mean price \eqn{MA_t} is \code{decay_factor} * \eqn{p_t}
#'     + \code{1-decay_factor} * \eqn{MA_{t-1}}.
#'
methods::setClass("LOAD", contains = "PortfolioSelector",
  slots = c(decay_factor = "numeric",
            regularization_factor = "numeric",
            time_window = "numeric",
            momentum_threshold = "numeric",
            wealth_factor_threshold = "numeric",
            prices = "numeric",
            price_means = "numeric")
)


#' @describeIn LOAD checks LOAD slots
#'
#' @title check_LOAD
#'
#' Ensures that the slots used for a \code{\linkS4class{LOAD}}
#' class are length 1 when they should be and that only the
#' initial prices variable is missing, if any
#'
#' @param object the \code{\linkS4class{LOAD}} instance to validate
#'
#' @return \code{TRUE} if validates the object, otherwise a character
#'     vector of the areas
#'
methods::setValidity("LOAD", function(object) {
  errors <- character()
  # make sure parameters have right length
  if(length(object@decay_factor) != 1) {
    msg <- glue::glue("decay_factor has length {length(object@decay_factor)}, \\
  should be length 1.")
    errors <- c(errors, msg)
  }
  if(length(object@regularization_factor) != 1) {
    msg <- glue::glue("regularization_factor has length \\
{length(object@regularization_factor)}, should be length 1.")
    errors <- c(errors, msg)
  }
  if(length(object@time_window) != 1) {
    msg <- glue::glue("time_window has length \\
{length(object@time_window)}, should be length 1.")
    errors <- c(errors, msg)
  }
  if(length(object@momentum_threshold) != 1) {
    msg <- glue::glue("momentum_threshold has length \\
{length(object@momentum_threshold)}, should be length 1.")
    errors <- c(errors, msg)
  }
  if(length(object@wealth_factor_threshold) != 1) {
    msg <- glue::glue("wealth_factor_threshold has length \\
{length(object@wealth_factor_threshold)}, should be length 1.")
    errors <- c(errors, msg)
  }
  if(any(is.na(c(object@decay_factor,
                 object@regularization_factor,
                 object@time_window,
                 object@momentum_threshold,
                 object@wealth_factor_threshold)))) {
    msg <- glue::glue("At least one of \\
(decay_factor, regularization_factor, time_window, momentum_threshold \\
wealth_factor_threshold) is missing.")
    errors <- c(errors, msg)
  }

  if(object@decay_factor < 0 || object@decay_factor > 1) {
    msg <- glue::glue("decay_factor of {object@decay_factor} not in [0,1].")
    errors <- c(errors, msg)
  }
  if(object@regularization_factor < 0) {
    msg <- glue::glue("regularization_factor must be non-negative.")
    errors <- c(errors, msg)
  }
  if(object@time_window < 2) {
    msg <- glue::glue("time_window must be at least 2")
    errors <- c(errors, msg)
  }
  if(!is_whole_number(object@time_window)) {
    msg <- glue::glue("time_window must be an integral value")
    errors <- c(errors, msg)
  }
  if(object@wealth_factor_threshold < 1) {
    msg <- glue::glue("wealth_factor_threshold must be at least 1.")
    errors <- c(errors, msg)
  }
  if(object@momentum_threshold < 0) {
    msg <- glue::glue("momentum_threshold must be non-negative.")
    errors <- c(errors, msg)
  }

  # make sure prices are the same size as price_relatives
  if(nrow(object@price_relatives) != nrow(object@prices)) {
    msg  <- glue::glue("price_relatives and prices have different number of rows")
    errors <- c(errors, msg)
  }
  if(ncol(object@price_relatives) != ncol(object@prices)) {
    msg  <- glue::glue("price_relatives and prices have different number of cols")
    errors <- c(errors, msg)
  }
  # make sure price_means are the right size
  if(nrow(object@price_relatives) != nrow(object@price_means)) {
    msg  <- glue::glue("price_relatives and price_means have different \\
number of rows")
    errors <- c(errors, msg)
  }
  if(ncol(object@price_relatives) != ncol(object@price_means)) {
    msg  <- glue::glue("price_relatives and price_means have different \\
number of cols")
    errors <- c(errors, msg)
  }

  if(length(errors) == 0) TRUE else errors
})


#' @describeIn LOAD initialize method for LOAD
#'
#' Ensures that if prices is not supplied, prices are computed
#' assuming initial prices are all 1.
#'
#' @importFrom magrittr %>%
#'
methods::setMethod("initialize", "LOAD",
function(.Object, decay_factor,
                  regularization_factor,
                  time_window,
                  momentum_threshold,
                  wealth_factor_threshold,
                  prices,
                  price_means,
                  ...) {
  .Object <- methods::callNextMethod(.Object, ...)
  # if prices is missing, assume initial prices are 1
  if(missing(prices)) {
    warn("prices is missing, this may lead to poor predictions.")
    nassets <- ncol(.Object@price_relatives)
    .Object@prices <- prices_from_relatives(.Object@price_relatives,
                                            rep(1, nassets))
  }
  if(missing(price_means)) {
    .Object@price_means <- .Object@prices %>%
      purrr::map(purrr::accumulate, ~ decay_factor *.y + (1-decay_factor) * .x)
  }
  methods::validObject(.Object)
})


#' @describeIn LOAD recommend the next portfolio using the LOAD system
#'
#' Recommends the next portfolio according to the LOAD system
#'
#' @param pf_selector The LOAD instance which is used to recommend a trade
#' @return The portfolio which should be traded to
#'
#' @importFrom magrittr %>%
#' @importFrom dplyr slice
#' @importFrom purrr map map_dbl imap_dbl
#' @importFrom glmnet glmnet
#'
methods::setMethod("recommend_next_portfolio", signature = "LOAD",
function(pf_selector) {
  # make sure we have the requisite number of trading periods
  if(pf_selector@trading_period < pf_selector@trading_window) {
    stop(glue::glue("LOAD requires at least trading_window-1 trading periods \\
before the current trading_period."))
  }
  tp <- pf_selector@trading_period
  # get predicted price based on regularized slope
  predicted_price(regularized_slope, asset) {
    if(regularized_slope > pf_selector@momentum_threshold) {
      pf_selector@prices %>%
        slice((tp - pf_selector@trading_window + 1):tp) %>%
        max() %>%
        return()
    }
    pf_selector@price_means[tp, asset]
  }
  # run ridge regression on each price
  # Note glmnet requires we have at least two variables but we only have
  # one, so we just duplicate the variable and sum them together
  predicted_price <- pf_selector@prices %>%
    slice((tp - pf_selector@trading_window + 1):tp) %>%
    map(glmnet,
        x = matrix(rep(1:pf_selector@trading_window, 2), ncol = 2),
        alpha = 0,  # alpha = 0 -> ridge regression
        lambda = pf_selector@regularization_factor) %>%
    map_dbl(`$`, beta) %>%
    map_dbl(sum) %>%
    imap_dbl(predicted_price)

  # use predicted price to get price relatives and compute
  # b_{t+1} according to the reference paper
  predicted_price_relatives <- predicted_price / pf_selector@prices[tp, ] %>%
    as.numeric()
  mean_zero_pred_pr <- predicted_price_relatives - mean(predicted_price_relatives)
  gamma <- as.numeric((
    pf_selector@wealth_factor_threshold
    - pf_selector@portfolio %*% predicted_price_relatives
  ) / (mean_zero_pred_pr %*% mean_zero_pred_pr))
  # if gamma <= 0 don't rebalance, otherwise do according to the
  # gamma * mean zero predicted price relatives
  if(gamma <= 0) {
    return(pf_selector@portfolio)
  }
  pf_selector@portfolio + gamma * mean_zero_pred_pr
})