用Rmosek求解Markowitz型投资组合优化问题

yrwegjxp  于 2023-09-27  发布在  其他
关注(0)|答案(1)|浏览(149)

我试图了解如何mosek工程和制定了一个制定了一个简单的最小xT * 西格玛 * x
我使用来自Rmosek码本的示例数据，并试图以圆锥曲线形式来表达问题。
按照mosek网站上的指南，二次项被t代替。然后使t最小化。

  1. # Example of input
  2. n      <- 8
  3. w      <- 1.0
  4. mu     <- c(0.07197349, 0.15518171, 0.17535435, 0.0898094, 0.42895777, 0.39291844, 0.32170722, 0.18378628)
  5. x0     <- c(0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)
  6. GT     <- rbind( c(0.30758, 0.12146, 0.11341, 0.11327, 0.17625, 0.11973, 0.10435, 0.10638),
  7.                  c(0.     , 0.25042, 0.09946, 0.09164, 0.06692, 0.08706, 0.09173, 0.08506),
  8.                  c(0.     , 0.     , 0.19914, 0.05867, 0.06453, 0.07367, 0.06468, 0.01914),
  9.                  c(0.     , 0.     , 0.     , 0.20876, 0.04933, 0.03651, 0.09381, 0.07742),
  10.                  c(0.     , 0.     , 0.     , 0.     , 0.36096, 0.12574, 0.10157, 0.0571 ),
  11.                  c(0.     , 0.     , 0.     , 0.     , 0.     , 0.21552, 0.05663, 0.06187),
  12.                  c(0.     , 0.     , 0.     , 0.     , 0.     , 0.     , 0.22514, 0.03327),
  13.                  c(0.     , 0.     , 0.     , 0.     , 0.     , 0.     , 0.     , 0.2202 ) )`
  15. library("Rmosek")
  17. # Sample covariance matrix
  18. Sigma <- t(GT)%*%GT
  20. # Number of assets
  21. n <- ncol(GT)
  23. # Problem definition
  24. prob <- list(sense="min")
  26. # Objective: We're minimizing t, which represents portfolio variance
  27. prob$c <- c(rep(0, n), 1)
  29. # Conic constraint: t >= x^T Sigma x
  30. prob$F <- rbind(
  31.   cbind(Matrix(0.0,nrow=1,ncol=n), 1), 
  32.   cbind(GT, Matrix(0, n, 1))
  33. )
  34. prob$g <- c(rep(0, n),1)
  36. prob$cones <- matrix(list("QUAD", 1+n, NULL), nrow=3, ncol=1)
  37. rownames(prob$cones) <- c("type","dim","conepar")
  39. # Linear constraint: Sum of x[i] = 1
  40. prob$A <- Matrix(c(rep(1, n), 0), 1)
  42. prob$bc <- rbind(blc=1,
  43.                  buc=1)
  45. # Bounds: 0 <= x[i] <= 1, 0 <= t
  46. prob$bx <- rbind(blx=c(rep(0, n + 1)), 
  47.                  bux=c(rep(1, n), Inf))
  49. # Solve the problem
  50. r <- mosek(prob,list(verbose=10))
  52. xTestMosekSimpleMinVarConic <- r$sol$itr$xx[1:n]
  54. stddevsol <- r$sol$itr$xx[n+1]
  56. stddevcalc <- (t(xTestMosekSimpleMinVarConic)%*%Sigma%*%xTestMosekSimpleMinVarConic)^0.5
  57. expret <- drop(mu %*% xTestMosekSimpleMinVarConic)
  59. print(sprintf("xTestMosekSimpleMinVarConic"))
  60. print(sprintf("Expected return: %.4e Std. deviation solver: %.4e Std. deviation calc: %.4e", round(expret,5), round(stddevsol,5), stddevcalc))
  61. print(round(xTestMosekSimpleMinVarConic,5))

当我解决这个问题时，结果是：

  1. [1] "Expected return: 1.7350e-01 Std. deviation solver: 1.0214e+00 Std. deviation calc: 2.0803e-01"
  2. [1] 0.13469 0.14807 0.29008 0.23863 0.00000 0.11060 0.07793 0.00000

然而，当我把同样的问题作为一个二次问题来解决时，我得到：

  1. [1] "Expected return: 1.6625e-01 Std. deviation: 2.0373e-01 "
  2. [1] 0.11332 0.11437 0.30221 0.18154 0.00000 0.05632 0.04519 0.18706

代码：

  1. library("Rmosek")
  3. Sigma <- t(GT)%*%GT
  4. # Number of assets
  5. n <- length(mu)
  7. # Define the QP problem
  8. prob <- list(sense = "min")
  10. # Convert the lower triangular part of the covariance matrix Sigma 
  11. # into the required list format for qobj
  12. i_indices <- c()
  13. j_indices <- c()
  14. values <- c()
  16. for (i in 1:n) {
  17.   for (j in 1:i) {
  18.     if (Sigma[i, j] != 0) {
  19.       i_indices <- append(i_indices, i)
  20.       j_indices <- append(j_indices, j)
  21.       values <- append(values, Sigma[i, j])
  22.     }
  23.   }
  24. }
  26. prob$qobj <- list(i = i_indices, j = j_indices, v = values)
  28. prob$c <-c(rep(0,n))
  30. # Equality constraint for sum of weights: 1^T * x = 1
  31. prob$A <- matrix(rep(1, n), 1)
  32. prob$bc <- rbind(blc = 1, buc = 1)
  34. # Non-negativity constraints on asset weights: x >= 0
  35. prob$bx <- rbind(blx = rep(0, n), bux = rep(Inf, n))
  37. # Solve the problem
  38. solution <- mosek(prob,list(verbose=1))
  40. xTestMosekSimpleQuad <- solution$sol$itr$xx
  41. stddev <- (t(xTestMosekSimpleQuad)%*%Sigma%*%xTestMosekSimpleQuad)^0.5
  42. expret <- drop(mu %*% xTestMosekSimpleQuad)
  44. print(sprintf("xTestMosekSimpleQuad"))
  45. print(sprintf("Expected return: %.4e Std. deviation: %.4e ", round(expret,5), round(stddev,5)))
  46. print(round(xTestMosekSimpleQuad,5))

当我使用标准求解器和Mosek求解器用CVXR解决问题时，得到的结果相同
CXVR代码：

  1. library("CVXR")
  3. Sigma <- t(GT)%*%GT
  5. ## Form problem
  6. w <- Variable(n)
  7. ret <- t(mu) %*% w
  8. risk <- quad_form(w, Sigma)
  9. constraints <- list(w >= 0, sum(w) == 1)
  10. objective <- risk
  11. prob <- Problem(Minimize(objective), constraints)
  12. result <- solve(prob, solver='MOSEK', verbose = TRUE)
  13. risk_data <- result$getValue(sqrt(risk))
  14. ret_data <- result$getValue(ret)
  15. w_data <- result$getValue(w)

CVXR将其作为圆锥问题处理：
问题名称：
客观意义：尽量减少
类型：圆锥优化问题（conic optimization problem）约束条件：19
仿射圆锥曲线：0
反意连词：0
锥体：1个
标量变量：19
基质变量：0
整型变量：0
优化程序已启动。已启动Presolve。线性相关性检查器已启动。线性相关性检查器已终止。清除器已启动。消除器中释放的约束：0消除器终止。消除器-尝试次数：1次：0.00
Lin. dep. -尝试：1次：0.00
Lin. dep. -原始尝试：1次成功：1
直线深度-双重尝试：0次成功：0
Lin. dep. -primal deps：0个双深度：0
预解决终止。时间：0.00
优化器-线程：12
优化器-已解决的问题：原始
优化器-约束：10
优化器-锥体：1
优化器-标量变量：18圆锥曲线：10
优化器-半定变量：0标量化：0
因子-设置时间：0.00
因子密集检测器。时间：0.00 GP订购时间：0.00
因子-因子前的非零值：54后因子：54
因素密集型调光：0次触发：1.20e +03
ITE PFEAS DFEAS GFEAS PRSTATUS POBJ DOBJ MU TIME
0 1.0e +00 6.0e-01 1.4e +00 0.00e +00 0.00000000e +00-3.973982284e-01 1.0e +00 0.00
1 4.1e-01 2.5e-01 1.1e-01 7.52e-01 4.724270646e-02-3.621647720e-01 4.1e-01 0.00
2 1.5e-01 8.9e-02 2.3e-02 2.20e +00 3.208365717e-02-4.069524762e-02 1.5e-01 0.00
3 2.0e-02 1.2e-02 1.1e-03 1.39e +00 4.349220651e-02 3.528925879e-02 2.0e-02 0.00
4 4.1e-03 2.4e-03 9.7e-05 1.05e +00 4.314371279e-02 4.153013503e-02 4.1e-03 0.00
5 6.2e-04 3.7e-04 7.0e-06 1.01e +00 4.148713867e-02 4.126481240e-02 6.2e-04 0.00
6 7.6e-05 4.6e-05 3.0e-07 1.00e +00 4.150471852e-02 4.147717911e-02 7.6e-05 0.00
7 1.0e-05 6.2e-06 1.5e-08 1.00e +00 4.150531100e-02 4.150163293e-02 1.0e-05 0.00
8 2.5e-07 1.5e-07 5.7e-11 1.00e +00 4.150536024e-02 4.150527034e-02 2.5e-07 0.00
9 5.1e-10 3.0e-10 5.0e-15 1.00e +00 4.150536306e-02 4.150536288e-02 5.0e-10 0.00
优化程序已终止。时间：0.00
内部点解决方案摘要问题状态：PRIMAL_AND_DUAL_FEASIBLE解决方案状态：最优原始。目标：4.1505363060e-02 nrm：1e +00 Viol. con：2e-10 var：0e +00个锥体：0e +00
双重目标：4.1505362883e-02 nrm：2e-01 Viol. con：0e +00 var：1e-10个锥体：0e +00
优化程序摘要优化程序-时间：0.00
内部点迭代：9次：0.00
基础识别-时间：0.00
原始-迭代：0时间：0.00
双迭代：0时间：0.00
清洁原始-迭代：0时间：0.00
干净的双重迭代：0时间：0.00
单纯形-时间：0.00
原始单纯形-迭代：0时间：0.00
对偶单纯形-迭代：0时间：0.00
混合整数-松弛：0时间：0.00
我在Mosek代码中的错误在哪里？

r

1条答案

按热度按时间
cuxqih21

cuxqih211#

我注意到的第一件事是你的评论说

  1. # Conic constraint: t >= x^T Sigma x

其中x^T Sigma x = x^T GT^T GT x = || GT x ||_2^2，但是你使用了二次锥，本质上得到

  1. t >= || GT x ||_2

其对应于t >= sqrt(x^T Sigma x)。如果t只出现在目标中，这种差异是好的（事实上，这是有益的，因为它为您提供了更好的缩放），但您需要在解释结果时意识到这种差异。
我注意到的第二件事是

  1. prob$g <- c(rep(0, n),1)

但这将常数1添加到圆锥约束的最后一个条目，本质上给出t >= || GT x + e_n ||_2，其中e_n是单位向量。尝试使用

  1. prob$g <- c(rep(0, n+1))

看看有没有用

展开查看全部
赞(0）分享  回复(0）举报  2023-09-27
我来回答

相关问题

    查看更多

    热门标签

    更多
    JavaquerypythonNode开发语言requestUtil数据库Table后端算法LoggerMessageElementParser

    最新问答

    更多
    • 蜀ICP备13028337号-1 大数据知识库 https://www.saoniuhuo.com © All rights reserved
    • 本站内容来源互联网,如果侵犯您的权益请联系我们删除, 联系方式：448109455@qq.com