GLPK/投资组合优化

本示例取自投资理论，属于金融经济学的分支。其目标是在已知历史表现的一组潜在投资中分配资本，以在收益和投资者的风险承受能力之间取得最佳平衡。

概要

本教程示例展示了如何使用 GLPK/GMPL 来优化投资组合，其中投资风险使用平均绝对偏差 (MAD) 标准来衡量。除了投资组合优化外，本示例还重点介绍了在其他 GMPL 应用中可能实用的技术，包括

Cholesky 分解正定矩阵。
从多元正态分布中生成样本，用于模拟历史收益。
从 GMPL 中为 GLPK 的随机数生成器设置种子。注意：新的randseed选项目前正在评估中。

相同作者的其他示例可在运筹学入门处找到。

平均绝对偏差 (MAD) 投资组合优化

投资组合优化指的是在金融资产集之间分配资本的过程，以实现收益和风险之间所需的平衡。该问题的经典马科维茨方法使用投资组合收益的预期方差来衡量风险。该标准为最佳投资组合中资产的相对权重提供了一个二次规划问题。

1991 年，Konno 和 Yamazaki ^[1] 提出了一个线性规划模型，用于投资组合优化，其中风险由平均绝对偏差 (MAD)（相对于预期收益的）来衡量。使用 MAD 作为风险度量标准产生的投资组合具有马科维茨投资组合不具备的几个理想特性，包括二阶随机支配。

如 Konno 和 Yamazaki 最初所提出的那样，首先需要一个收益历史记录 $R_{i}(t_{n})$ ，用于资产集中的每项资产 $S$ 。时间 $t_{n}$ 的收益由资产价格的变化决定， $R_{i}(t_{n})={(P_{i}(t_{n})-P_{i}(t_{n-1}))}/{P_{i}(t_{n-1})}$ 。对于每项资产，预期收益通过以下公式估算：

{\bar {R}}_{i}\approx {\frac {1}{N}}\sum _{n=1}^{N}R_{i}(t_{n})

.

投资者指定一个最低要求收益 $R_{p}$ 。投资组合优化问题是确定分配给每项资产的总投资的比例，即 $w_{i}$ ，其最小化与平均值的平均绝对偏差

\min _{w_{i}}{\frac {1}{N}}\sum _{n=1}^{N}|\sum _{i\in S}w_{i}(R_{i}(t_{n})-{\bar {R}}_{i})|

受所需回报和固定总投资约束

{\begin{array}{rcl}\sum _{i\in S}w_{i}{\bar {R}}_{i}&\geq &R_{p}\\\quad \\\sum _{i\in S}w_{i}&=&1\end{array}}

最小所需回报的值， $R_{p}$ 表达了投资者的风险承受能力。 $R_{p}$ 的值越小，则可行解空间越大，导致投资组合的 MAD 风险指标值越低。增加 $R_{p}$ 会导致投资组合具有更高的风险。风险和回报之间的关系是投资理论的基本原理。

此公式不会对权重 $w_{i}$ 设置个体边界。特别是， $w_{i}<0$ 对应于相关资产的卖空。如果投资者对卖空、单一资产的最大投资额或投资组合多元化的其他要求设置限制，则可以添加约束条件。根据可用投资机会的集合，额外的约束条件可能会导致不可行的投资问题。

实现细节

本示例的目的是展示在 GMPL 中对投资组合优化问题的建模。主要任务是将平均绝对偏差目标转化为线性形式，并模拟给定一组潜在投资的均值和协方差的的历史回报。

MAD 目标的重新表述

目标函数的以下公式改编自 Gerald Curnuejols 和 Reha Tutuncu (2007) ^[2] 的 "金融优化方法"，该公式反过来遵循 Feinstein 和 Thapa (1993) ^[3]。通过引入决策变量 $y_{n}\geq 0$ 和 $z_{n}\geq 0$ ， $n=1,\ldots ,N$ ，简化模型使目标变为

\min _{w_{i},y_{n},z_{n}}{\frac {1}{N}}\sum _{n=1}^{N}(y_{n}+z_{n})

受制于约束

{\begin{array}{rcl}\sum _{i\in S}w_{i}{\bar {R}}_{i}&\geq &R_{p}\\\\\sum _{i\in S}w_{i}&=&1\\\\y_{n}-z_{n}&=&\sum _{i\in S}w_{i}(R_{i}(t_{n})-{\bar {R}}_{i})\quad {\mbox{for}}\quad n=1,\ldots ,N\end{array}}

正如 Feinstein 和 Thapa 所讨论的，这个版本将问题简化为 $N+2$ 个约束，在 $2N+{\mbox{card}}(S)$ 个决策变量中，注意到 card 函数。

为 GLPK 伪随机数生成器设置种子

不幸的是，GMPL 没有提供一种方法来设置 GLPK 中的伪随机数生成器的种子。相反，以下 GMPL 代码片段使用 gmtime() 函数来查找自 1970 年以来的秒数。丢弃前导数字可以避免后续的溢出错误，而平方可以返回一个每秒变化很大的数字。提取最低位数字可以产生一个介于 0 到 100,000 之间的数字，该数字决定在后续使用之前对 GLPK 的伪随机数生成器进行采样的次数。这种技巧没有保证其统计特性。

 /* A workaround for the lack of a way to seed the PRNG in GMPL */
param utc := prod {1..2} (gmtime()-1000000000);
param seed := utc - 100000*floor(utc/100000);
check sum{1..seed} Uniform01() > 0;

历史回报的模拟

对于这个实现，历史回报是通过假设了解资产回报的均值和协方差来模拟的。我们从一个平均回报向量开始 ${\bar {R}}$ 和协方差矩阵 $\Sigma$ ，由下式估计

\Sigma _{ij}\approx {\frac {1}{N-1}}\sum _{n=1}^{N}(R_{i}(t_{n})-{\bar {R}}_{i})(R_{j}(t_{n})-{\bar {R}}_{j})

历史回报的模拟需要从多变量正态分布中生成样本。为此，我们计算 Cholesky 分解，其中，对于一个正半定 $\Sigma$ ， $\Sigma =CC^{T}$ ，并且 $C$ 是一个下三角矩阵。以下 GMPL 代码片段实现了 Cholesky-Crout 算法。

/* Cholesky Lower Triangular Decomposition of the Covariance Matrix */
param C{i in S, j in S : i >= j} :=
    if i = j then
        sqrt(Sigma[i,i]-(sum {k in S : k < i} (C[i,k]*C[i,k])))
    else
        (Sigma[i,j]-sum{k in S : k < j} C[i,k]*C[j,k])/C[j,j];

如果没有错误检查，除非 $\Sigma$ 是正定的，否则该代码片段将失败。对于真实世界的数据，协方差矩阵通常是正定的，所以这通常不是问题。但是，如果尝试在投资资产集中包含无风险资产（如政府债券），这就会成为一个问题。

一旦计算出 Cholesky 因子 $C$ ，模拟收益向量 $R(t_{n})$ 由 $R(t_{n})={\bar {R}}+CZ(t_{n})$ 给出，其中 $Z(t_{n})$ 的元素是从均值为零、方差为一的正态分布中独立抽取的样本。

/* Simulated returns */
param N default 5000;
set T := 1..N;
param R{i in S, t in T} := Rbar[i] + sum {j in S : j <= i} C[i,j]*Normal(0,1);

PortfolioMAD.mod

将此模型保存为PortfolioMAD.mod.

/* Portfolio Optimization using Mean Absolute Deviation

Jeff Kantor 
December 4, 2009
Revised: December 6, 2009 to fix problem with random variate generation.
Revised: December 7, 2009 to add a 'seeding' of the PRNG
Revised: July 8, 2010 reformatted for GLPK Wikibook */

/* Stock Data */

set S;                                              # Set of stocks
param Rbar{S};                                      # Means of projected returns
param Sigma{S,S};                                   # Covariance of projected returns
param Rp  default (1/card(S))*sum{i in S} Rbar[i];  # Lower bound on portfolio return

/* Generate sample data */

/* Cholesky Lower Triangular Decomposition of the Covariance Matrix */
param C{i in S, j in S : i >= j} :=
    if i = j then
        sqrt(Sigma[i,i]-(sum {k in S : k < i} (C[i,k]*C[i,k])))
    else
        (Sigma[i,j]-sum{k in S : k < j} C[i,k]*C[j,k])/C[j,j];

 /* A workaround for the lack of a way to seed the PRNG in GMPL */
param utc := prod {1..2} (gmtime()-1000000000);
param seed := utc - 100000*floor(utc/100000);
check sum{1..seed} Uniform01() > 0;

/* Simulated returns */
param N default 5000;
set T := 1..N;
param R{i in S, t in T} := Rbar[i] + sum {j in S : j <= i} C[i,j]*Normal(0,1);

/* MAD Optimization */

var w{S};                     # Portfolio Weights with Bounds
var y{T} >= 0;                # Positive deviations (non-negative)
var z{T} >= 0;                # Negative deviations (non-negative)

minimize MAD: (1/card(T))*sum {t in T} (y[t] + z[t]);

s.t. C1: sum {s in S} w[s]*Rbar[s] >= Rp;
s.t. C2: sum {s in S} w[s] = 1;
s.t. C3 {t in T}: (y[t] - z[t]) = sum{s in S} (R[s,t]-Rbar[s])*w[s];

solve;

/* Report */

/* Input Data */
printf "\n\nStock Data\n\n";
printf "         Return   Variance\n";
printf {i in S} "%5s   %7.4f   %7.4f\n", i, Rbar[i], Sigma[i,i];

printf "\nCovariance Matrix\n\n";
printf "     ";
printf {j in S} " %7s ", j;
printf "\n";
for {i in S} {
    printf "%5s  " ,i;
    printf {j in S} " %7.4f ", Sigma[i,j];
    printf "\n";
}

/* MAD Optimal Portfolio */
printf "\n\nMinimum Absolute Deviation (MAD) Portfolio\n\n";
printf "  Return   = %7.4f\n",Rp;
printf "  Variance = %7.4f\n\n", sum {i in S, j in S} w[i]*w[j]*Sigma[i,j];
printf "         Weight\n";
printf {s in S} "%5s   %7.4f\n", s, w[s];
printf "\n";

data;

/* Data for monthly returns on four selected stocks for a three
year period ending December 4, 2009 */

# Simulation Horizon
param N := 5000;

# Minimum acceptable investment return
param Rp := 0.01;

# Historical returns on assets
param : S : Rbar :=
    AAPL    0.0308
    GE     -0.0120
    GS      0.0027
    XOM     0.0018 ;
    
# Covariance on asset returns
param   Sigma :
            AAPL    GE      GS      XOM  :=
    AAPL    0.0158  0.0062  0.0088  0.0022
    GE      0.0062  0.0136  0.0064  0.0011
    GS      0.0088  0.0064  0.0135  0.0008
    XOM     0.0022  0.0011  0.0008  0.0022 ;

end;

结果

典型的输出如下。结果表明，与单个资产的风险回报表现相比，精心设计的投资组合可以表现出明显改善的风险回报表现。

Stock Data

         Return   Variance
 AAPL    0.0308    0.0158
   GE   -0.0120    0.0136
   GS    0.0027    0.0135
  XOM    0.0018    0.0022

Covariance Matrix

         AAPL       GE       GS      XOM 
 AAPL    0.0158   0.0062   0.0088   0.0022
   GE    0.0062   0.0136   0.0064   0.0011
   GS    0.0088   0.0064   0.0135   0.0008
  XOM    0.0022   0.0011   0.0008   0.0022

Minimum Absolute Deviation (MAD) Portfolio

  Return   =  0.0100
  Variance =  0.0036

         Weight
 AAPL    0.2794
   GE    0.0002
   GS    0.1120
  XOM    0.6084

可能的扩展

此示例有许多扩展，可以使其在现实世界中更加有用

添加 Cholesky 分解的错误检查。
添加权重系数的上限和下限（例如，没有卖空）约束。
添加为投资组合指定无风险资产的能力。
添加约束以实施多元化要求。
添加风险回报权衡的参数分析（需要外部脚本）。
开发一种从 GMPL 内部对 GLPK 的随机数生成器进行播种的更好方法。
添加一个ODBC 数据库接口，用于经验（而不是模拟）历史收益数据。

参考

↑ Konno, Hiroshi; Yamazaki, Hiroaki (1991). "Mean-absolute deviation portfolio optimization model and its applications to Tokyo stock market". Management Science. 37: 519–531.
↑ Curnuejols, Gerald; Tutuncu, Reha (2007). Optimization Methods in Finance. Cambridge University Press.
↑ Feinstein, Charles D.; Thapa, Mukund N. (1993). "A Reformulation of a Mean-Absolute Deviation Portfolio Optimization Model". Management Science. 39: 1552–1553.

[1] Konno, Hiroshi; Yamazaki, Hiroaki (1991). "Mean-absolute deviation portfolio optimization model and its applications to Tokyo stock market". Management Science. 37: 519–531.

[2] Curnuejols, Gerald; Tutuncu, Reha (2007). Optimization Methods in Finance. Cambridge University Press.

[3] Feinstein, Charles D.; Thapa, Mukund N. (1993). "A Reformulation of a Mean-Absolute Deviation Portfolio Optimization Model". Management Science. 39: 1552–1553.

[1]

[2]

[3]