使用“R”中的“fPortfolio”包找到给定风险水平下具有最大回报的投资组合

Question

我正在尝试使用

fPortfolio

中的

包找到在给定风险水平下具有最大回报的投资组合。然而，该软件包中似乎存在一个错误，即在更改目标风险水平时，最佳投资组合（基于最大回报）不会改变，这表明它没有考虑用户指定的目标风险水平。其他人已发现此问题（请参阅下面的链接），我在下面提供了一个最小的可重现示例。据我所知，这个问题尚未在package文档（存档于here）或作者关于fPortfolio的书中（存档于here）中得到解决。该错误已在近 10 年前被描述过，但仍未得到修复，因此我认为我们可以向社区求助以帮助解决该问题。我创建此问题的目标是总结已知内容，识别包代码的潜在问题，并希望获得帮助来解决它。

以下是关于此问题的 StackOverflow 问题：

在 R 包 fPortfolio 中设置目标风险（已存档此处）
在 R 包 fPortfolio 中使用目标风险或目标回报（已存档此处）

这是一个最小的可重现示例：

library("quantmod")
#> Warning: package 'quantmod' was built under R version 4.3.3
#> Loading required package: xts
#> Loading required package: zoo
#> 
#> Attaching package: 'zoo'
#> The following objects are masked from 'package:base':
#> 
#>     as.Date, as.Date.numeric
#> Loading required package: TTR
#> Warning: package 'TTR' was built under R version 4.3.3
#> Registered S3 method overwritten by 'quantmod':
#>   method            from
#>   as.zoo.data.frame zoo
library("fPortfolio")
#> Warning: package 'fPortfolio' was built under R version 4.3.3
#> Loading required package: timeDate
#> Warning: package 'timeDate' was built under R version 4.3.2
#> Loading required package: timeSeries
#> Warning: package 'timeSeries' was built under R version 4.3.3
#> 
#> Attaching package: 'timeSeries'
#> The following object is masked from 'package:zoo':
#> 
#>     time<-
#> The following objects are masked from 'package:graphics':
#> 
#>     lines, points
#> Loading required package: fBasics
#> Warning: package 'fBasics' was built under R version 4.3.3
#> 
#> Attaching package: 'fBasics'
#> The following object is masked from 'package:TTR':
#> 
#>     volatility
#> Loading required package: fAssets
#> Warning: package 'fAssets' was built under R version 4.3.3
library("dplyr")
#> Warning: package 'dplyr' was built under R version 4.3.2
#> 
#> ######################### Warning from 'xts' package ##########################
#> #                                                                             #
#> # The dplyr lag() function breaks how base R's lag() function is supposed to  #
#> # work, which breaks lag(my_xts). Calls to lag(my_xts) that you type or       #
#> # source() into this session won't work correctly.                            #
#> #                                                                             #
#> # Use stats::lag() to make sure you're not using dplyr::lag(), or you can add #
#> # conflictRules('dplyr', exclude = 'lag') to your .Rprofile to stop           #
#> # dplyr from breaking base R's lag() function.                                #
#> #                                                                             #
#> # Code in packages is not affected. It's protected by R's namespace mechanism #
#> # Set `options(xts.warn_dplyr_breaks_lag = FALSE)` to suppress this warning.  #
#> #                                                                             #
#> ###############################################################################
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:timeSeries':
#> 
#>     filter, lag
#> The following objects are masked from 'package:xts':
#> 
#>     first, last
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union

# Download historical stock prices
symbols <- c("AAPL", "MSFT", "GOOGL", "AMZN", "JNJ", "JPM", "V", "PG", "XOM", "TSLA")
quantmod::getSymbols(symbols)
#>  [1] "AAPL"  "MSFT"  "GOOGL" "AMZN"  "JNJ"   "JPM"   "V"     "PG"    "XOM"  
#> [10] "TSLA"

# Calculate stock returns
prices <- do.call(merge, lapply(symbols, function(sym) quantmod::Cl(get(sym))))
returns <- na.omit(TTR::ROC(prices, type = "discrete"))

returns_ts <- timeSeries::as.timeSeries(returns)

# Create portfolio
portfolioSpec <- fPortfolio::portfolioSpec()

# Identify tangency portfolio
tangencyPortfolio <- fPortfolio::tangencyPortfolio(
  data = returns_ts,
  spec = portfolioSpec)

# Extract optimal weights for tangency portfolio
fPortfolio::getWeights(tangencyPortfolio)
#>  AAPL.Close  MSFT.Close GOOGL.Close  AMZN.Close   JNJ.Close   JPM.Close 
#>  0.21136483  0.12495822  0.01223435  0.18743248  0.00000000  0.00000000 
#>     V.Close    PG.Close   XOM.Close  TSLA.Close 
#>  0.23495738  0.06078310  0.00000000  0.16826965

# Define target risk levels
targetRisks <- seq(0, 0.3, by = 0.01)

# Initialize storage for optimal portfolios
optimalPortfolios <- list()
optimalWeights_list <- list()

# Find optimal weightings for each target risk level
for (risk in targetRisks) {
  # Create a portfolio optimization specification with the target risk
  portfolioSpec <- fPortfolio::portfolioSpec()
  fPortfolio::setTargetRisk(portfolioSpec) <- risk
  
  # Solve for the maximum return at this target risk
  optimal_portfolio <- fPortfolio::maxreturnPortfolio(
    returns_ts,
    spec = portfolioSpec)
  
  # Store the optimal portfolio
  optimalPortfolios[[as.character(risk)]] <- optimal_portfolio
  
  # Store the optimal portfolio weights with risk level
  optimal_weights <- fPortfolio::getWeights(optimal_portfolio)
  optimalWeights_list[[as.character(risk)]] <- c(RiskLevel = risk, optimal_weights)
}

optimalWeightsByRisk <- dplyr::bind_rows(optimalWeights_list)
optimalWeightsByRisk
#> # A tibble: 31 × 11
#>    RiskLevel AAPL.Close MSFT.Close GOOGL.Close AMZN.Close JNJ.Close JPM.Close
#>        <dbl>      <dbl>      <dbl>       <dbl>      <dbl>     <dbl>     <dbl>
#>  1      0        0.0276          0      0.0282     0.0466     0.393         0
#>  2      0.01     0.0276          0      0.0282     0.0466     0.393         0
#>  3      0.02     0.0276          0      0.0282     0.0466     0.393         0
#>  4      0.03     0.0276          0      0.0282     0.0466     0.393         0
#>  5      0.04     0.0276          0      0.0282     0.0466     0.393         0
#>  6      0.05     0.0276          0      0.0282     0.0466     0.393         0
#>  7      0.06     0.0276          0      0.0282     0.0466     0.393         0
#>  8      0.07     0.0276          0      0.0282     0.0466     0.393         0
#>  9      0.08     0.0276          0      0.0282     0.0466     0.393         0
#> 10      0.09     0.0276          0      0.0282     0.0466     0.393         0
#> # ℹ 21 more rows
#> # ℹ 4 more variables: V.Close <dbl>, PG.Close <dbl>, XOM.Close <dbl>,
#> #   TSLA.Close <dbl>

sessionInfo()
#> R version 4.3.1 (2023-06-16 ucrt)
#> Platform: x86_64-w64-mingw32/x64 (64-bit)
#> Running under: Windows 11 x64 (build 22631)
#> 
#> Matrix products: default
#> 
#> 
#> locale:
#> [1] LC_COLLATE=English_United States.utf8 
#> [2] LC_CTYPE=English_United States.utf8   
#> [3] LC_MONETARY=English_United States.utf8
#> [4] LC_NUMERIC=C                          
#> [5] LC_TIME=English_United States.utf8    
#> 
#> time zone: America/Chicago
#> tzcode source: internal
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#>  [1] dplyr_1.1.4         fPortfolio_4023.84  fAssets_4023.85    
#>  [4] fBasics_4032.96     timeSeries_4032.109 timeDate_4032.109  
#>  [7] quantmod_0.4.26     TTR_0.24.4          xts_0.14.0         
#> [10] zoo_1.8-12         
#> 
#> loaded via a namespace (and not attached):
#>  [1] utf8_1.2.4          generics_0.1.3      bitops_1.0-7       
#>  [4] robustbase_0.99-3   slam_0.1-51         lattice_0.22-6     
#>  [7] digest_0.6.36       magrittr_2.0.3      rneos_0.4-0        
#> [10] evaluate_0.24.0     grid_4.3.1          fCopulae_4022.85   
#> [13] mvtnorm_1.2-5       fastmap_1.2.0       jsonlite_1.8.8     
#> [16] energy_1.7-11       fansi_1.0.6         kernlab_0.9-32     
#> [19] truncnorm_1.0-9     XML_3.99-0.17       numDeriv_2016.8-1.1
#> [22] ecodist_2.1.3       mnormt_2.1.1        cli_3.6.3          
#> [25] rlang_1.1.4         Rglpk_0.6-5.1       gsl_2.1-8          
#> [28] reprex_2.1.1        withr_3.0.0         yaml_2.3.9         
#> [31] tools_4.3.1         parallel_4.3.1      mvnormtest_0.1-9-3 
#> [34] boot_1.3-30         curl_5.2.1          R6_2.5.1           
#> [37] vctrs_0.6.5         stats4_4.3.1        lifecycle_1.0.4    
#> [40] fs_1.6.4            MASS_7.3-60.0.1     Rsolnp_1.16        
#> [43] pkgconfig_2.0.3     pillar_1.9.0        fMultivar_4031.84  
#> [46] glue_1.7.0          Rcpp_1.0.13         tidyselect_1.2.1   
#> [49] tibble_3.2.1        DEoptimR_1.1-3      xfun_0.46          
#> [52] rstudioapi_0.16.0   knitr_1.48          spatial_7.3-17     
#> [55] htmltools_0.5.8.1   igraph_2.0.3        rmarkdown_2.27     
#> [58] compiler_4.3.1      quadprog_1.5-8      sn_2.1.1           
#> [61] RCurl_1.98-1.16

^{创建于 2024-07-23，使用 reprex v2.1.1}

这里是

maxreturnPortfolio()

函数的封装代码：

> fPortfolio::maxreturnPortfolio
function (data, spec = portfolioSpec(), constraints = "LongOnly") 
{
    data = portfolioData(data, spec)
    if (is.null(getTargetRisk(spec))) {
        stop("Missing target risk for maximum return optimization.")
    }
    else {
        Solver = match.fun(getSolver(spec))
        portfolio = Solver(data, spec, constraints)
        setWeights(spec) = portfolio$weights
        setStatus(spec) = portfolio$status
        Title = "Return Maximized Efficient Portfolio"
    }
    portfolio <- feasiblePortfolio(data, spec, constraints)
    portfolio@call <- match.call()
    portfolio@title <- Title
    portfolio
}

它似乎调用了

feasiblePortfolio()

函数：

> fPortfolio::feasiblePortfolio
function (data, spec = portfolioSpec(), constraints = "LongOnly") 
{
    Data <- portfolioData(data, spec)
    if (inherits(data, "fPFOLIODATA")) 
        data <- getSeries(Data)
    assetsNames <- getUnits(Data)
    Spec <- spec
    Constraints <- portfolioConstraints(Data, spec, constraints)
    if (is.null(getWeights(spec))) {
        stop("Missing weights")
    }
    weights <- as.vector(getWeights(spec))
    names(weights) <- assetsNames
    if (inherits(getSeries(Data), "timeSeries")) {
        targetReturn <- c(mean = (Data@statistics$mean %*% weights)[[1]], 
            mu = (Data@statistics$mu %*% weights)[[1]])
        setTargetReturn(spec) <- targetReturn
        Cov <- Data@statistics$Cov
        cov <- sqrt((weights %*% Cov %*% weights)[[1]])
        if (getType(spec) == "SPS") {
            myCheck <- TRUE
            funSigma <- match.fun(getObjective(spec)[1])
            rcov <- funSigma(as.vector(weights))
        }
        else {
            Sigma <- Data@statistics$Sigma
            rcov <- sqrt((weights %*% Sigma %*% weights)[[1]])
        }
        alpha <- getAlpha(spec)
        returns <- getDataPart(getSeries(Data)) %*% weights
        VaR <- quantile(returns, alpha, type = 1)
        CVaR <- VaR - 0.5 * mean(((VaR - returns) + abs(VaR - 
            returns)))/alpha
        targetRisk <- c(cov, rcov, -CVaR, -VaR)
        names(targetRisk) <- c("Cov", "Sigma", "CVaR", "VaR")
        alpha <- getAlpha(Spec)
    }
    else if (inherits(getSeries(Data), "logical")) {
        targetReturn <- c(mean = (Data@statistics$mean %*% weights)[[1]], 
            mu = NA)
        setTargetReturn(spec) <- targetReturn
        Cov <- Data@statistics$Cov
        cov <- sqrt((weights %*% Cov %*% weights)[[1]])
        targetRisk <- c(cov, NA, NA, NA)
        names(targetRisk) <- c("Cov", "Sigma", "CVaR", "VaR")
        alpha <- NA
    }
    covRiskBudgets <- (weights * Cov %*% weights)[, 1]/cov^2
    names(covRiskBudgets) <- assetsNames
    Portfolio <- new("fPFOLIOVAL", portfolio = list(weights = weights, 
        covRiskBudgets = covRiskBudgets, targetReturn = targetReturn, 
        targetRisk = targetRisk, targetAlpha = alpha, status = getStatus(spec)))
    new("fPORTFOLIO", call = match.call(), data = Data, spec = Spec, 
        constraints = Constraints, portfolio = Portfolio, title = "Feasible Portfolio", 
        description = description())
}

据我所知，根据代码（我不是专家），

feasiblePortfolio()

函数似乎创建并分配目标风险，而不是接受它作为用户的输入（即，

 targetRisk <- c(cov, rcov, -CVaR, -VaR)

），虽然我不确定。正如一位 SO 用户指出的那样，“[它] 似乎与

efficientPortfolio()

函数和

spec

中使用的求解器有关。显然，无论您输入什么

targetRisk

，求解器都会输入一个与您打印投资组合时得到的

Cov

相对应的目标”。

任何有助于识别和解决此问题的帮助将不胜感激。

Answer 1

如here所述，当允许卖空时，使用

target risk

优化

maxreturnPortfolio

似乎可以按预期工作，定义

setSolver(portfolioSpec)= "solveRshortExact"

。对于

"LongOnly"

优化，正如您所说，问题会持续很长时间，在这种情况下

maxreturnPortfolio

总是返回最小无功投资组合优化作为解决方案，这是错误的。

我针对这种情况提出了一种基于

portfolioFrontier

解决方案的解决方法。该想法考虑前沿风险中最接近目标风险的两个点，并对两个投资组合进行线性组合以迭代至所需的

tolerance

风险。

# Function to aproximated portfolio for a target risk (tar_risk)
# based on efficiente frontier estimation (Eff_front) by linear combination
# of two nearest risks in the frontier, for most convex problems
# like efficient frontier will work within a 'tolerance',
# for more restricted problems a verification of further restrictions
# must be implemented
# for "LongOnly" restrictions

Aprox_maxreturn_Portfolio <- function(tar_risk,Eff_front,tolerance){
  # get parameters from Eff_front
  covarest = getCov(Eff_front) 
  Expret = getMean(Eff_front)
  sigma= Eff_front@portfolio@portfolio$targetRisk[,"Sigma"]
  fret= Eff_front@portfolio@portfolio$targetReturn[,"mean"]
  fweights = Eff_front@portfolio@portfolio$weights
  n=ncol(fweights) # assets
  EFpoints = fret >= fret[which.min(sigma)] # the efficient points only
  fweights = fweights[EFpoints,]
  sigma = sigma[EFpoints]
  ransigma = range(sigma)
  if (!(tar_risk>=ransigma[1] & tar_risk<= ransigma[2]))
    stop("tar_risk is not inside efficient frontier risk interval ", round(ransigma,6),"\n")
  posrisk= which(diff(findInterval(sigma, tar_risk))==1)
  # iterate whithin a tolerance
  iniws = fweights
  diffrisk = 1e9
  while( diffrisk > tolerance){
    # linear combin adjacent ef portfolios
    wz = seq(0,1,length.out=100)
    fe=matrix(rep(0,n*length(wz)),nrow = length(wz), ncol=n, byrow=TRUE)
    colnames(fe)=colnames(fweights)
    for (i in 1:n) fe[,i]=wz*iniws[posrisk,i]+(1-wz)*iniws[posrisk+1,i]
    ExpRetPort=rep(0,length(wz)) #
    RiskPort=ExpRetPort
    for (j in 1:length(wz)){ ExpRetPort[j]=fe[j,]%*%Expret # exp ret
    specj = portfolioSpec(); setTargetReturn(specj) = ExpRetPort[j]
    effportj =  efficientPortfolio(getData(Eff_front), spec = specj, getConstraints(Eff_front))
    RiskPort[j] = getTargetRisk(effportj)["Sigma"]
    #RiskPort[j]=(t(fe[j,])%*%covarest%*%fe[j,])^.5 # risk
    } 
    posrisk= which(abs(diff(findInterval(RiskPort, tar_risk)))==1)
    diffrisk= abs(diff(RiskPort[c(posrisk,posrisk+1)]))
    iniws = fe
  }
  
  # compute final portfolio
  wopt=(RiskPort[c(posrisk,posrisk+1)]-tar_risk)/abs(diff(RiskPort[c(posrisk,posrisk+1)]))
  wopt = abs(wopt[1])
  tport =wopt*fe[posrisk,]+(1-wopt)*fe[posrisk+1,]
  tret=tport%*%Expret # exp ret
  tar_riskc=(t(tport)%*%covarest%*%tport)^.5 
  
  # solution
  list(portfolio=tport,return = tret, risk=tar_riskc )
}

示例：

# Create portfolio
portfolioSpec <- portfolioSpec()

Constraints = "LongOnly"
setNFrontierPoints(portfolioSpec) = 80 
Data = returns_ts

## Eff Frontier
Eff_front = portfolioFrontier(Data,portfolioSpec, Constraints)

#plot
frontierPlot(Eff_front, pch = 19, title =FALSE)
Eff_front@portfolio@portfolio$minriskPortfolio

c(min(Eff_front@portfolio@portfolio$targetRisk[,"Sigma"]),
  max(Eff_front@portfolio@portfolio$targetRisk[,"Sigma"]))
# solve a problem
tar_risk = 0.0091
tolerance = 1e-7

Aprox_maxreturn_Portfolio(tar_risk,Eff_front,tolerance)

$portfolio
 AAPL.Close  MSFT.Close GOOGL.Close  AMZN.Close   JNJ.Close   JPM.Close     V.Close    PG.Close   XOM.Close 
 0.04379226  0.00000000  0.02848604  0.05911637  0.36892972  0.00000000  0.02692962  0.36391607  0.09255214 
 TSLA.Close 
 0.01627778 

$return
             [,1]
[1,] 0.0004797394

$risk
            [,1]
[1,] 0.009100002

多种目标风险模拟

# simulate for several points
# Define target risk levels

targetRisks <- seq(0.009064718,0.035888956, by = 0.002)

# Initialize storage for optimal portfolios
optimalPortfolios <- list()
optimalWeights_list <- list()

# Find optimal weightings for each target risk level

for (risk in targetRisks) {
  # Create a portfolio optimization specification with the target risk
  portfolio1 = Aprox_maxreturn_Portfolio(risk,Eff_front,tolerance)
  # Store the optimal portfolio
  optimalPortfolios[[as.character(risk)]] <- portfolio1
  # Store the optimal portfolio weights with risk level
  optimal_weights <- portfolio1$portfolio
  optimalWeights_list[[as.character(risk)]] <- c(targetRisk = portfolio1$risk,
                                                 targetReturn=portfolio1$return,
                                                 targetRiskd  = risk,
                                                 optimal_weights)
}
optimalWeightsByRisk <- dplyr::bind_rows(optimalWeights_list)
optimalWeightsByRisk

#MAE
mean(abs(apply(optimalWeightsByRisk[,c(1,3)],1,diff)))

#[1] 2.612458e-08

绘制两个边界

# plot
upperFrontier <- frontierPoints(Eff_front, frontier = "upper", 
                                return = "mean", risk = "Sigma", auto = T)
dataefffront = rbind(data.frame(upperFrontier,Frontier="Efficient"),
                     data.frame(optimalWeightsByRisk[,1:2],Frontier="Approx"))
library(ggplot2)
ggplot(dataefffront,aes(x=targetRisk,y=targetReturn,colour = Frontier))+
  geom_point()

使用“R”中的“fPortfolio”包找到给定风险水平下具有最大回报的投资组合

问题描述投票：0回答：1

1个回答

最新问题

使用“R”中的“fPortfolio”包找到给定风险水平下具有最大回报的投资组合

问题描述 投票：0回答：1

1个回答

最新问题

问题描述投票：0回答：1