R оптимизация купить продать

Я считаю, что это можно смоделировать как небольшую модель смешанного целочисленного программирования (MIP).

Вот реализация с использованием CVXR:

> library(CVXR)
> 
> # data
> price = c(12, 11, 12, 13, 16, 17, 18, 17, 18, 16, 17, 13)
> capacity = 25
> max_units_buy = 4
> max_units_sell = 8
> 
> # number of time periods
> NT <- length(price)
> 
> # Decision variables
> inv = Variable(NT,integer=T)
> buy = Variable(NT,integer=T)
> sell = Variable(NT,integer=T)
> 
> # Lag operator
> L = cbind(rbind(0,diag(NT-1)),0)
> 
> # optimization model
> problem <- Problem(Maximize(sum(price*(sell-buy))),
+                    list(inv == L %*% inv + buy - sell,
+                         inv >= 0, inv <= capacity,
+                         buy >= 0, buy <= max_units_buy,
+                         sell >= 0, sell <= max_units_sell))
> result <- solve(problem,verbose=T)
GLPK Simplex Optimizer, v4.47
84 rows, 36 columns, 119 non-zeros
*     0: obj =  0.000000000e+000  infeas = 0.000e+000 (12)
*    35: obj = -1.040000000e+002  infeas = 0.000e+000 (0)
OPTIMAL SOLUTION FOUND
GLPK Integer Optimizer, v4.47
84 rows, 36 columns, 119 non-zeros
36 integer variables, none of which are binary
Integer optimization begins...
+    35: mip =     not found yet >=              -inf        (1; 0)
+    35: >>>>> -1.040000000e+002 >= -1.040000000e+002   0.0% (1; 0)
+    35: mip = -1.040000000e+002 >=     tree is empty   0.0% (0; 1)
INTEGER OPTIMAL SOLUTION FOUND
> cat("status:",result$status)
status: optimal
> cat("objective:",result$value)
objective: 104
> print(result$getValue(buy))
      [,1]
 [1,]    4
 [2,]    4
 [3,]    4
 [4,]    4
 [5,]    4
 [6,]    0
 [7,]    0
 [8,]    4
 [9,]    0
[10,]    4
[11,]    0
[12,]    0
> print(result$getValue(sell))
      [,1]
 [1,]    0
 [2,]    0
 [3,]    0
 [4,]    0
 [5,]    0
 [6,]    8
 [7,]    8
 [8,]    0
 [9,]    8
[10,]    0
[11,]    4
[12,]    0
> print(result$getValue(inv))
      [,1]
 [1,]    4
 [2,]    8
 [3,]   12
 [4,]   16
 [5,]   20
 [6,]   12
 [7,]    4
 [8,]    8
 [9,]    0
[10,]    4
[11,]    0
[12,]    0
> 

Добавлена ​​возможность иметь начальный инвентарь и создана функция для поэтапной оптимизации с учетом покупки / продажи, зависящей от уровня запасов.

library(tidyverse)
library(lubridate)
library(CVXR)

init_fce <- function(.df_storage, .df_bounds, .type = "max"){
  if(.type == "max"){
    .df_storage$max_buy <- max(.df_bounds$max_buy)
    .df_storage$max_sell <- max(.df_bounds$max_sell)
  } else if(.type == "min"){
    .df_storage$max_buy <- min(.df_bounds$max_buy)
    .df_storage$max_sell <- min(.df_bounds$max_sell)
  } else if(.type == "mean"){
    .df_storage$max_buy <- mean(.df_bounds$max_buy)
    .df_storage$max_sell <- mean(.df_bounds$max_sell)
  }

  .df_storage
}
optim_fce <- function(.df){
  
  # Decision variables
  m_inv_tot = Variable(nrow(.df), integer = T)
  m_buy = Variable(nrow(.df), integer = T)
  m_sell = Variable(nrow(.df), integer = T)
  # Lag operator
  m_L = cbind(rbind(0, diag(nrow(.df) - 1)), 0)
  
  objetive <- Maximize(sum(.df$price*(m_sell-m_buy)))
  
  constraints <- list(
    m_inv_tot == m_L %*% m_inv_tot + .df$inv_init + m_buy - m_sell, # L %*% result$getValue(inv) + result$getValue(buy) - result$getValue(sell)
    m_inv_tot >= 0, m_inv_tot <= .df$capacity,
    m_buy >= 0, m_buy <= .df$max_buy,
    m_sell >= 0, m_sell <= .df$max_sell
  )
  
  problem <- Problem(objetive, constraints)
  result <- solve(problem) # , verbose=T
  
  .df <- .df %>%
    mutate(
      buy = (result$getValue(m_buy) %>% as.vector()),
      sell = (result$getValue(m_sell)  %>% as.vector()),
      inventory_real = (result$getValue(m_inv_tot) %>% as.vector())
    )
  
  .df
}
set_limits_fce <- function(.df_storage, .df_bounds){
  .df_storage <- .df_storage %>%
    select(-max_buy, -max_sell) %>%
    mutate(capacity_usage_pct_prec = lag(inventory_real, default = inv_init[1])/capacity) %>%
    crossing(.df_bounds %>% select(-segment)) %>%
    filter(capacity_usage_pct_prec >= lbound, capacity_usage_pct_prec < ubound) %>%
    mutate(
      within_bounds = (buy <= max_buy) & (sell <= max_sell) 
    ) %>%
    select(-lbound, -ubound)

  .df_storage
}
get_results <- function(.df_storage){
  if( any(!.df_storage$within_bounds) ){
    print("result not within bounds")
  } else{
    .df_storage$profit <- .df_storage$sell * .df_storage$price - .df_storage$buy * .df_storage$price
    print(sum(.df_storage$profit))
  }
  
  .df_storage
}


A1_storage <- tibble(
  date = ymd("2020-06-01") + months(0:11),
  price = c(12, 11, 12, 13, 16, 17, 18, 17, 18, 16, 17, 13),
  inv_init = c(3, rep(0, 11)),
  capacity = 25
)

A2_bounds <- tibble(
  segment = c("0%-30%", "30%-65%", "65%-70%", "70%-100%"),
  lbound = c(0, 0.3, 0.65, 0.7), 
  ubound = c(0.3, 0.65, 0.7, 1),
  max_buy = c(4,3,2,2),
  max_sell = c(4,6,6,8)
)

B1_max <- init_fce(A1_storage, A2_bounds, .type = "max") %>%
  optim_fce() %>%
  set_limits_fce(.df_bounds = A2_bounds) %>%
  get_results() %>%
  optim_fce() %>%
  set_limits_fce(.df_bounds = A2_bounds) %>%
  get_results() %>%
  optim_fce() %>%
  set_limits_fce(.df_bounds = A2_bounds) %>%
  get_results() %>%
  optim_fce() %>%
  set_limits_fce(.df_bounds = A2_bounds) %>%
  get_results()