R中的递归optim（）函数导致错误
Recursive optim() function in R causes errors
我正在尝试使用 R 中的 optim() 函数来最小化矩阵运算的值。在这种情况下，我试图将一组股票的波动性降至最低，这些股票的个体 returns 相互协变。被最小化的 objective 函数是 calculate_portfolio_variance.
library(quantmod)
filter_and_sort_symbols <- function(symbols)
{
# Name: filter_and_sort_symbols
# Purpose: Convert to uppercase if not
# and remove any non valid symbols
# Input: symbols = vector of stock tickers
# Output: filtered_symbols = filtered symbols
# convert symbols to uppercase
symbols <- toupper(symbols)
# Validate the symbol names
valid <- regexpr("^[A-Z]{2,4}$", symbols)
# Return only the valid ones
return(sort(symbols[valid == 1]))
}
# Create the list of stock tickers and check that they are valid symbols
tickers <- filter_and_sort_symbols(c("AAPL", "NVDA", "MLM", "AA"))
benchmark <- "SPY"
# Set the start and end dates
start_date <- "2007-01-01"
end_date <- "2019-01-01"
# Gather the stock data using quantmod library
getSymbols(Symbols=tickers, from=start_date, to=end_date, auto.assign = TRUE)
getSymbols(benchmark, from=start_date, to=end_date, auto.assign = TRUE)
# Create a matrix of only the adj. prices
price_matrix <- NULL
for(ticker in tickers){price_matrix <- cbind(price_matrix, get(ticker)[,6])}
# Set the column names for the price matrix
colnames(price_matrix) <- tickers
benchmark_price_matrix <- NULL
benchmark_price_matrix <- cbind(benchmark_price_matrix, get(benchmark)[,6])
# Compute log returns
returns_matrix <- NULL
for(ticker in tickers){returns_matrix <- cbind(returns_matrix, annualReturn(get(ticker)))}
returns_covar <- cov(returns_matrix)
colnames(returns_covar) <- tickers
rownames(returns_covar) <- tickers
# get average returns for tickers and benchmark
ticker_avg <- NULL
for(ticker in tickers){ticker_avg <- cbind(ticker_avg, colMeans(annualReturn(get(ticker))))}
colnames(ticker_avg) <- tickers
benchmark_avg <- colMeans(annualReturn(get(benchmark)))
# create the objective function
calculate_portfolio_variance <- function(allocations, returns_covar, ticker_avg, benchmark_avg)
{
# Name: calculate_portfolio_variance
# Purpose: Computes expected portfolio variance, to be used as the minimization objective function
# Input: allocations = vector of allocations to be adjusted for optimality; returns_covar = covariance matrix of stock returns
# ticker_avg = vector of average returns for all tickers, benchmark_avg = benchmark avg. return
# Output: Expected portfolio variance
# get benchmark volatility
benchmark_variance <- (sd(annualReturn(get(benchmark))))^2
# scale allocations for 100% investment
allocations <- as.matrix(allocations/sum(allocations))
# get the naive allocations
naive_allocations <- rep(c(1/ncol(ticker_avg)), times=ncol(ticker_avg))
portfolio_return <- sum(t(allocations)*ticker_avg)
portfolio_variance <- t(allocations)%*%returns_covar%*%allocations
# constraints = portfolio expected return must be greater than benchmark avg. return and
# portfolio variance must be less than benchmark variance (i.e. a better reward at less risk)
if(portfolio_return < benchmark_avg | portfolio_variance > benchmark_variance)
{
allocations <- naive_allocations
}
portfolio_variance <- t(allocations)%*%returns_covar%*%allocations
return(portfolio_variance)
}
# Specify lower and upper bounds for the allocation percentages
lower <- rep(0, ncol(returns_matrix))
upper <- rep(1, ncol(returns_matrix))
# Initialize the allocations by evenly distributing among all tickers
set.seed(1234)
allocations <- rep(1/length(tickers), times=length(tickers))
当我手动调用 objective 函数时，它 returns 一个预期的值：
> calculate_portfolio_variance(allocations, returns_covar, ticker_avg, benchmark_avg)
[,1]
[1,] 0.1713439
但是，当我使用 optim() 函数时 returns 出现错误：
> optim_result <- optim(par=allocations, fn=calculate_portfolio_variance(allocations, ticker_avg, benchmark_avg), lower=lower, upper=upper, method="L-BFGS-B")
Error in t(allocations) %*% returns_covar : non-conformable arguments
我不确定原因，但可能与 optim() 递归使用 allocations 变量的方式有关。我该怎么做才能解决这个问题？
编辑：FWIW，其他优化策略有效（差分进化、模拟退火）但我更喜欢使用梯度下降，因为它要快得多
如果将第一个参数重命名为 par 并且您将应用 t() 的顺序切换到该侧翼矩阵乘法运算中使用的参数向量，则不会发生错误：
cpv <- function(par, returns_covar=returns_covar, ticker_avg, benchmark_avg)
{
# Name: calculate_portfolio_variance
# Purpose: Computes expected portfolio variance, to be used as the minimization objective function
# Input: allocations = vector of allocations to be adjusted for optimality; returns_covar = covariance matrix of stock returns
# ticker_avg = vector of average returns for all tickers, benchmark_avg = benchmark avg. return
# Output: Expected portfolio variance
# get benchmark volatility
benchmark_variance <- (sd(annualReturn(get(benchmark))))^2
# scale allocations for 100% investment
par <- as.matrix(par/sum(par))
# get the naive allocations
naive_allocations <- rep(c(1/ncol(ticker_avg)), times=ncol(ticker_avg))
portfolio_return <- sum(t(par)*ticker_avg);print(par)
portfolio_variance <- t(par)%*%returns_covar%*%par
# constraints = portfolio expected return must be greater than benchmark avg. return and
# portfolio variance must be less than benchmark variance (i.e. a better reward at less risk)
if(portfolio_return < benchmark_avg | portfolio_variance > benchmark_variance)
{
par <- naive_allocations
}
portfolio_variance <- t(par)%*%returns_covar%*%par
return(portfolio_variance)
}
我在代码中留下了par的调试打印，显示运行它的结果的顶部
optim_result <- optim(par=allocations, fn=cpv, lower=lower, upper=upper, returns_covar=returns_covar, ticker_avg=ticker_avg, benchmark_avg=benchmark_avg, method="L-BFGS-B")
[,1]
[1,] 0.25
[2,] 0.25
[3,] 0.25
[4,] 0.25
[,1]
[1,] 0.2507493
[2,] 0.2497502
[3,] 0.2497502
[4,] 0.2497502
[,1]
[1,] 0.2492492
[2,] 0.2502503
[3,] 0.2502503
[4,] 0.2502503
#--- snipped output of six more iterations.
...结果：
> optim_result
$par
[1] 0.25 0.25 0.25 0.25
$value
[1] 0.1713439
$counts
function gradient
1 1
$convergence
[1] 0
$message
[1] "CONVERGENCE: NORM OF PROJECTED GRADIENT <= PGTOL"
正如我在对一个无关问题的评论中所说，optim 函数首先尝试升高然后降低 par 中的第一个元素，然后尝试对第二个、第三个和第四个元素执行相同的操作。在那一点上没有发现任何改进，它“决定”它收敛于局部最小值并宣布收敛。
我应该指出 optim 的代码是 rather old and the author of the original algorithm, Dr Nash, has placed an updated version on CRAN in the form of the optimx package。他说 optim 当时很好，但他认为如果不成功，应该尝试其他程序。
我正在尝试使用 R 中的 optim() 函数来最小化矩阵运算的值。在这种情况下，我试图将一组股票的波动性降至最低，这些股票的个体 returns 相互协变。被最小化的 objective 函数是 calculate_portfolio_variance.
library(quantmod)
filter_and_sort_symbols <- function(symbols)
{
# Name: filter_and_sort_symbols
# Purpose: Convert to uppercase if not
# and remove any non valid symbols
# Input: symbols = vector of stock tickers
# Output: filtered_symbols = filtered symbols
# convert symbols to uppercase
symbols <- toupper(symbols)
# Validate the symbol names
valid <- regexpr("^[A-Z]{2,4}$", symbols)
# Return only the valid ones
return(sort(symbols[valid == 1]))
}
# Create the list of stock tickers and check that they are valid symbols
tickers <- filter_and_sort_symbols(c("AAPL", "NVDA", "MLM", "AA"))
benchmark <- "SPY"
# Set the start and end dates
start_date <- "2007-01-01"
end_date <- "2019-01-01"
# Gather the stock data using quantmod library
getSymbols(Symbols=tickers, from=start_date, to=end_date, auto.assign = TRUE)
getSymbols(benchmark, from=start_date, to=end_date, auto.assign = TRUE)
# Create a matrix of only the adj. prices
price_matrix <- NULL
for(ticker in tickers){price_matrix <- cbind(price_matrix, get(ticker)[,6])}
# Set the column names for the price matrix
colnames(price_matrix) <- tickers
benchmark_price_matrix <- NULL
benchmark_price_matrix <- cbind(benchmark_price_matrix, get(benchmark)[,6])
# Compute log returns
returns_matrix <- NULL
for(ticker in tickers){returns_matrix <- cbind(returns_matrix, annualReturn(get(ticker)))}
returns_covar <- cov(returns_matrix)
colnames(returns_covar) <- tickers
rownames(returns_covar) <- tickers
# get average returns for tickers and benchmark
ticker_avg <- NULL
for(ticker in tickers){ticker_avg <- cbind(ticker_avg, colMeans(annualReturn(get(ticker))))}
colnames(ticker_avg) <- tickers
benchmark_avg <- colMeans(annualReturn(get(benchmark)))
# create the objective function
calculate_portfolio_variance <- function(allocations, returns_covar, ticker_avg, benchmark_avg)
{
# Name: calculate_portfolio_variance
# Purpose: Computes expected portfolio variance, to be used as the minimization objective function
# Input: allocations = vector of allocations to be adjusted for optimality; returns_covar = covariance matrix of stock returns
# ticker_avg = vector of average returns for all tickers, benchmark_avg = benchmark avg. return
# Output: Expected portfolio variance
# get benchmark volatility
benchmark_variance <- (sd(annualReturn(get(benchmark))))^2
# scale allocations for 100% investment
allocations <- as.matrix(allocations/sum(allocations))
# get the naive allocations
naive_allocations <- rep(c(1/ncol(ticker_avg)), times=ncol(ticker_avg))
portfolio_return <- sum(t(allocations)*ticker_avg)
portfolio_variance <- t(allocations)%*%returns_covar%*%allocations
# constraints = portfolio expected return must be greater than benchmark avg. return and
# portfolio variance must be less than benchmark variance (i.e. a better reward at less risk)
if(portfolio_return < benchmark_avg | portfolio_variance > benchmark_variance)
{
allocations <- naive_allocations
}
portfolio_variance <- t(allocations)%*%returns_covar%*%allocations
return(portfolio_variance)
}
# Specify lower and upper bounds for the allocation percentages
lower <- rep(0, ncol(returns_matrix))
upper <- rep(1, ncol(returns_matrix))
# Initialize the allocations by evenly distributing among all tickers
set.seed(1234)
allocations <- rep(1/length(tickers), times=length(tickers))
当我手动调用 objective 函数时，它 returns 一个预期的值：
> calculate_portfolio_variance(allocations, returns_covar, ticker_avg, benchmark_avg)
[,1]
[1,] 0.1713439
但是，当我使用 optim() 函数时 returns 出现错误：
> optim_result <- optim(par=allocations, fn=calculate_portfolio_variance(allocations, ticker_avg, benchmark_avg), lower=lower, upper=upper, method="L-BFGS-B")
Error in t(allocations) %*% returns_covar : non-conformable arguments
我不确定原因，但可能与 optim() 递归使用 allocations 变量的方式有关。我该怎么做才能解决这个问题？
编辑：FWIW，其他优化策略有效（差分进化、模拟退火）但我更喜欢使用梯度下降，因为它要快得多
如果将第一个参数重命名为 par 并且您将应用 t() 的顺序切换到该侧翼矩阵乘法运算中使用的参数向量，则不会发生错误：
cpv <- function(par, returns_covar=returns_covar, ticker_avg, benchmark_avg)
{
# Name: calculate_portfolio_variance
# Purpose: Computes expected portfolio variance, to be used as the minimization objective function
# Input: allocations = vector of allocations to be adjusted for optimality; returns_covar = covariance matrix of stock returns
# ticker_avg = vector of average returns for all tickers, benchmark_avg = benchmark avg. return
# Output: Expected portfolio variance
# get benchmark volatility
benchmark_variance <- (sd(annualReturn(get(benchmark))))^2
# scale allocations for 100% investment
par <- as.matrix(par/sum(par))
# get the naive allocations
naive_allocations <- rep(c(1/ncol(ticker_avg)), times=ncol(ticker_avg))
portfolio_return <- sum(t(par)*ticker_avg);print(par)
portfolio_variance <- t(par)%*%returns_covar%*%par
# constraints = portfolio expected return must be greater than benchmark avg. return and
# portfolio variance must be less than benchmark variance (i.e. a better reward at less risk)
if(portfolio_return < benchmark_avg | portfolio_variance > benchmark_variance)
{
par <- naive_allocations
}
portfolio_variance <- t(par)%*%returns_covar%*%par
return(portfolio_variance)
}
我在代码中留下了par的调试打印，显示运行它的结果的顶部
optim_result <- optim(par=allocations, fn=cpv, lower=lower, upper=upper, returns_covar=returns_covar, ticker_avg=ticker_avg, benchmark_avg=benchmark_avg, method="L-BFGS-B")
[,1]
[1,] 0.25
[2,] 0.25
[3,] 0.25
[4,] 0.25
[,1]
[1,] 0.2507493
[2,] 0.2497502
[3,] 0.2497502
[4,] 0.2497502
[,1]
[1,] 0.2492492
[2,] 0.2502503
[3,] 0.2502503
[4,] 0.2502503
#--- snipped output of six more iterations.
...结果：
> optim_result
$par
[1] 0.25 0.25 0.25 0.25
$value
[1] 0.1713439
$counts
function gradient
1 1
$convergence
[1] 0
$message
[1] "CONVERGENCE: NORM OF PROJECTED GRADIENT <= PGTOL"
正如我在对一个无关问题的评论中所说，optim 函数首先尝试升高然后降低 par 中的第一个元素，然后尝试对第二个、第三个和第四个元素执行相同的操作。在那一点上没有发现任何改进，它“决定”它收敛于局部最小值并宣布收敛。
我应该指出 optim 的代码是 rather old and the author of the original algorithm, Dr Nash, has placed an updated version on CRAN in the form of the optimx package。他说 optim 当时很好，但他认为如果不成功，应该尝试其他程序。