空売り制約下での分散最小化問題を解く、ということで、今回のテーマは最適化問題。

最適化問題を解く上で使用するライブラリcvxpyについてまとめる。

前回の平均分散アプローチによる最小分散フロンティアを求めた時との違いは、空売り制約がかかったことで、解析解がなくなったこと。解析解がない→最適化問題として解くということなのかと。

やりたいこと

コード

import numpy as np
import numpy.linalg as lin
import cvxpy as cvx

Mu=np.array([1.0,3.0,1.5,6.0,4.5])#期待値
Std=np.array([5.0,10.0,7.5,15.0,11.0]) #分散
CorrMatrix = np.array([[1.00, 0.25, 0.18, 0.10, 0.25],
                       [0.25, 1.00, 0.36, 0.20, 0.20],
                       [0.18, 0.36, 1.00, 0.25, 0.36],
                       [0.10, 0.20, 0.25, 1.00, 0.45],
                       [0.25, 0.20, 0.36, 0.45, 1.00]])

Sigma=np.diag(Std).dot(CorrMatrix).dot(np.diag(Std))#分散共分散行列
iota=np.ones(Mu.shape)
inv_Sigma=lin.inv(Sigma)

#空売り制約下での分散最小化問題の設定
Weight=cvx.Variable(Mu.shape[0])
Target_Return=cvx.Parameter(nonneg=True)
Risk_Variance=cvx.quad_form(Weight,Sigma)
constraints=[Weight.T*Mu == Target_Return,
             cvx.sum(Weight)==1.0,
             Weight>=0.0]
Opt_Portfolio=cvx.Problem(cvx.Minimize(Risk_Variance),constraints)

15行目までは、前回分散最小ポートフォリオを求めたものと同じ。

17行目から、cvxpyを用いた最適化問題を設定する

cvxpy

最適化問題の設定と、それを解く二段階でメモ。

最適化問題の設定

最適化問題なので、まずは目的関数と制約条件の記述が必須。

その２つを表しているのが、21行目の関数cvx.Problem

cvx.Problem('目的関数',[制約条件])といった感じで表現。

目的関数は、Sigma(分散共分散行列)とWeight(各資産の比率)であらわされるので、まだ定義されていない変数Weightを定義。また、各資産比率は、この最適化問題の中で制御する変数であるので、cvxpyの関数Variableで定義。それが18行目の以下。

Weight=cvx.Variable(Mu.shape[0])

列数を指定しないとき、変数は列ベクトルとなる。21行目コード中の制約条件の中で、Weightの転置行列とMuの積はWeight.T*Muと表現されているのが気になっていたがこれで納得した。

前回の.dotを用いた表現だと、Weight.dot(Mu)といったように、左から乗じたWeightは断りなく行ベクトルとして扱われていたが、今回はvariable関数によって列ベクトルとされているので、転置が必要ということではないかと。

また、今回の最適化問題では、パラメータとしてTarget_Return(期待収益率)を設定。

それを表しているのが、19行目の以下。

Target_Return=cvx.Parameter(nonneg=True)

オプションのnonnegはパラメータの符号は正という意味。

そして、目的関数を20行目で定義。

quad_formは二次形式を計算するcvxpy関数。二次形式については略。こんなのやった覚えない。。。

いざ解く

まずは準備。

いざ解く上では、まずパラメータを設定する必要あり。以下の配列を作成。

V_Target=np.linspace(Mu.min(),Mu.max(),num=250)

あとは、それぞれのパラメータに対して最適化問題を解いた結果の制御変数Weightと分散Risk_Varianceを格納するための配列を以下で作成。

V_Risk=np.zeros(V_Target.shape)
V_Weight=np.zeros((V_Target.shape[0],Mu.shape[0]))

いよいよ解く。

V_Target=np.linspace(Mu.min(),Mu.max(),num=250)
V_Risk=np.zeros(V_Target.shape)
V_Weight=np.zeros((V_Target.shape[0],Mu.shape[0]))

for idx,Target_Return.value in enumerate(V_Target):
    Opt_Portfolio.solve()
    V_Weight[idx,:]=Weight.value.T
    V_Risk[idx]=np.sqrt(Risk_Variance.value)

cvxpyでの関数で作成した変数の値を取り出すときは.valueを使うみたい

cvxpyまとめ

最適化問題設定

最適化問題なので、以下の２つは必須。

制御変数の設定cvx.Variable
最適化オブジェクトの設定cvx.Problem

必要に応じて

パラメータの設定cvx.Parameter

解く

最適化オブジェクト.solve()

最後に

今回cvxpyを用いて作成したコード

import numpy as np
import numpy.linalg as lin
import cvxpy as cvx
import matplotlib.pyplot as plt

Mu=np.array([1.0,3.0,1.5,6.0,4.5])#期待値
Std=np.array([5.0,10.0,7.5,15.0,11.0]) #分散
CorrMatrix = np.array([[1.00, 0.25, 0.18, 0.10, 0.25],
                       [0.25, 1.00, 0.36, 0.20, 0.20],
                       [0.18, 0.36, 1.00, 0.25, 0.36],
                       [0.10, 0.20, 0.25, 1.00, 0.45],
                       [0.25, 0.20, 0.36, 0.45, 1.00]])

Sigma=np.diag(Std).dot(CorrMatrix).dot(np.diag(Std))#分散共分散行列
iota=np.ones(Mu.shape)
inv_Sigma=lin.inv(Sigma)

#空売り制約下での分散最小化問題の設定
Weight=cvx.Variable(Mu.shape[0])
Target_Return=cvx.Parameter(nonneg=True)
Risk_Variance=cvx.quad_form(Weight,Sigma)
constraints=[Weight.T*Mu == Target_Return,
             cvx.sum(Weight)==1.0,
             Weight>=0.0]
Opt_Portfolio=cvx.Problem(cvx.Minimize(Risk_Variance),constraints)

V_Target=np.linspace(Mu.min(),Mu.max(),num=250)
V_Risk=np.zeros(V_Target.shape)
V_Weight=np.zeros((V_Target.shape[0],Mu.shape[0]))

for idx,Target_Return.value in enumerate(V_Target):
    Opt_Portfolio.solve()
    V_Weight[idx,:]=Weight.value.T
    V_Risk[idx]=np.sqrt(Risk_Variance.value)

#可視化
fig1=plt.figure(1,facecolor='w')
plt.plot(V_Risk,V_Target,'k-')
plt.xlabel(u'vol[%]')
plt.ylabel(u'Return')

plt.show()

結果

やっとできて嬉しい、、

ぶぶぶろぐ

ライブラリcvxpy