![](zzz/editor/img/code.gif)
# -*- coding: utf-8 -*-
# NumPyの読み込み
import numpy as np
# Pandasの読み込み
import pandas as pd
# Cvxoptの読み込み
import cvxopt as co
# MatplotlibのPyplotモジュールの読み込み
import matplotlib.pyplot as plt
# 日本語フォントの設定
# 最小分散ポートフォリオの計算
def MV_Portfolio(ReturnData, TargetReturn):
## ReturnData: 収益率データ
## TargetReturn: 目標収益率
## Output: 最適ポートフォリオ
co.solvers.options['show_progress'] = False
MeanVector = np.mean(ReturnData, axis=0)
Deviation = ReturnData - MeanVector
NumberOfPeriods = ReturnData.shape[0]
NumberOfAssets = ReturnData.shape[1]
NumberOfElements = NumberOfAssets + NumberOfPeriods
P = co.spmatrix(1.0/NumberOfPeriods, \
range(NumberOfAssets, NumberOfElements), \
range(NumberOfAssets, NumberOfElements))
q = co.matrix(0.0, (NumberOfElements, 1))
G = co.spmatrix(-1.0, range(NumberOfAssets), range(NumberOfAssets), \
(NumberOfAssets, NumberOfElements))
h = co.matrix(0.0, (NumberOfAssets, 1))
A = co.sparse([ \
co.sparse([ \
[co.matrix(Deviation)], \
[co.spmatrix(-1.0, range(NumberOfPeriods), range(NumberOfPeriods))]\
]), \
co.sparse([ \
[co.matrix(np.c_[MeanVector, np.ones(NumberOfAssets)]).trans()], \
[co.spmatrix([], [], [], (2, NumberOfPeriods))] \
]) \
])
b = co.matrix(np.r_[np.zeros(NumberOfPeriods), TargetReturn, 1.0])
return co.solvers.qp(P, q, G, h, A, b)
# 最小絶対偏差ポートフォリオの計算
def AD_Portfolio(ReturnData, TargetReturn):
## ReturnData: 収益率データ
## TargetReturn: 目標収益率
## Output: 最適ポートフォリオ
co.solvers.options['show_progress'] = False
MeanVector = np.mean(ReturnData, axis=0)
Deviation = ReturnData - MeanVector
NumberOfPeriods = ReturnData.shape[0]
NumberOfAssets = ReturnData.shape[1]
NumberOfElements = NumberOfAssets + NumberOfPeriods
c = co.matrix([co.matrix(0.0, (NumberOfAssets, 1)), \
co.matrix(1.0/NumberOfPeriods, (NumberOfPeriods, 1))])
G = co.sparse([ \
co.spmatrix(-1.0, range(NumberOfPeriods), \
range(NumberOfAssets, NumberOfElements)), \
co.sparse([ \
[co.matrix(-Deviation)], \
[co.spmatrix(-1.0, range(NumberOfPeriods), range(NumberOfPeriods))]\
]), \
co.spmatrix(-1.0, range(NumberOfAssets), range(NumberOfAssets), \
(NumberOfAssets, NumberOfElements)) \
])
h = co.matrix(0.0, (2*NumberOfPeriods + NumberOfAssets, 1))
A = co.sparse([ \
[co.matrix(np.c_[MeanVector, np.ones(NumberOfAssets)]).trans()], \
[co.spmatrix([], [], [], (2, NumberOfPeriods))] \
])
b = co.matrix([TargetReturn, 1.0])
return co.solvers.lp(c, G, h, A, b)
# 最適期待ショートフォール・ポートフォリオの計算
def ES_Portfolio(ReturnData, TargetReturn, Alpha):
## ReturnData: 収益率データ
## TargetReturn: 目標収益率
## Alpha: 下方確率
## Output: 最適ポートフォリオ
co.solvers.options['show_progress'] = False
MeanVector = np.mean(ReturnData, axis=0)
NumberOfPeriods = ReturnData.shape[0]
NumberOfAssets = ReturnData.shape[1]
NumberOfElements = NumberOfAssets + NumberOfPeriods
c = co.matrix( \
[co.matrix(0.0, (NumberOfAssets, 1)), \
co.matrix(1.0/(Alpha*NumberOfPeriods), (NumberOfPeriods, 1)), -1.0])
G = co.sparse([ \
co.spmatrix(-1.0, range(NumberOfPeriods), \
range(NumberOfAssets, NumberOfElements), \
(NumberOfPeriods, NumberOfElements + 1)), \
co.sparse([ \
[co.matrix(-ReturnData)], \
[co.spmatrix(-1.0, range(NumberOfPeriods), range(NumberOfPeriods))], \
[co.matrix(1.0, (NumberOfPeriods, 1))] \
]), \
co.spmatrix(-1.0, range(NumberOfAssets), range(NumberOfAssets), \
(NumberOfAssets, NumberOfElements + 1)) \
])
h = co.matrix(0.0, (2*NumberOfPeriods + NumberOfAssets, 1))
A = co.sparse([ \
[co.matrix(np.c_[MeanVector, np.ones(NumberOfAssets)]).trans()], \
[co.spmatrix([], [], [], (2, NumberOfPeriods + 1))] \
])
b = co.matrix([TargetReturn, 1.0])
return co.solvers.lp(c, G, h, A, b)
# 株価データの読み込み
stockvalue = pd.read_csv('stocks.csv', index_col=0)
# 最小分散フロンティアの計算
R = (stockvalue.diff()/stockvalue.shift(1))[1:] * 100.0
T = R.shape[0]
mu = R.mean().values
Sigma = (R.cov().values)*((T-1.0)/T)
V_Target = np.linspace(mu.min(), mu.max(), num=50)
V_SD = np.array( \
[np.sqrt(2.0*MV_Portfolio(R.values, Target)['primal objective']) \
for Target in V_Target])
V_AD = np.array( \
[2.0*AD_Portfolio(R.values, Target)['primal objective'] \
for Target in V_Target])
V_ES = np.array( \
[ES_Portfolio(R.values, Target, 0.1)['primal objective'] \
for Target in V_Target])
# フロンティアのグラフの作成
()
fig1 = plt.figure(1, facecolor='w')
plt.subplot(1, 3, 1)
plt.plot(V_SD, V_Target, 'b-')
plt.plot(np.sqrt(np.diagonal(Sigma)), mu, 'rx')
plt.xlabel(u'standard deviation')
plt.ylabel(u'Expected rate of return')
plt.subplot(1, 3, 2)
plt.plot(V_AD, V_Target, 'b-')
plt.plot((R - R.mean()).abs().mean().values, mu, 'rx')
plt.xlabel(u'Absolute deviation')
plt.subplot(1, 3, 3)
plt.plot(V_ES, V_Target, 'b-')
plt.plot((-R[R <= R.quantile(0.1)]).mean().values, mu, 'rx')
plt.xlabel(u'Expected shortfall')
# NumPyの読み込み
import numpy as np
# Pandasの読み込み
import pandas as pd
# Cvxoptの読み込み
import cvxopt as co
# MatplotlibのPyplotモジュールの読み込み
import matplotlib.pyplot as plt
# 日本語フォントの設定
# 最小分散ポートフォリオの計算
def MV_Portfolio(ReturnData, TargetReturn):
## ReturnData: 収益率データ
## TargetReturn: 目標収益率
## Output: 最適ポートフォリオ
co.solvers.options['show_progress'] = False
MeanVector = np.mean(ReturnData, axis=0)
Deviation = ReturnData - MeanVector
NumberOfPeriods = ReturnData.shape[0]
NumberOfAssets = ReturnData.shape[1]
NumberOfElements = NumberOfAssets + NumberOfPeriods
P = co.spmatrix(1.0/NumberOfPeriods, \
range(NumberOfAssets, NumberOfElements), \
range(NumberOfAssets, NumberOfElements))
q = co.matrix(0.0, (NumberOfElements, 1))
G = co.spmatrix(-1.0, range(NumberOfAssets), range(NumberOfAssets), \
(NumberOfAssets, NumberOfElements))
h = co.matrix(0.0, (NumberOfAssets, 1))
A = co.sparse([ \
co.sparse([ \
[co.matrix(Deviation)], \
[co.spmatrix(-1.0, range(NumberOfPeriods), range(NumberOfPeriods))]\
]), \
co.sparse([ \
[co.matrix(np.c_[MeanVector, np.ones(NumberOfAssets)]).trans()], \
[co.spmatrix([], [], [], (2, NumberOfPeriods))] \
]) \
])
b = co.matrix(np.r_[np.zeros(NumberOfPeriods), TargetReturn, 1.0])
return co.solvers.qp(P, q, G, h, A, b)
# 最小絶対偏差ポートフォリオの計算
def AD_Portfolio(ReturnData, TargetReturn):
## ReturnData: 収益率データ
## TargetReturn: 目標収益率
## Output: 最適ポートフォリオ
co.solvers.options['show_progress'] = False
MeanVector = np.mean(ReturnData, axis=0)
Deviation = ReturnData - MeanVector
NumberOfPeriods = ReturnData.shape[0]
NumberOfAssets = ReturnData.shape[1]
NumberOfElements = NumberOfAssets + NumberOfPeriods
c = co.matrix([co.matrix(0.0, (NumberOfAssets, 1)), \
co.matrix(1.0/NumberOfPeriods, (NumberOfPeriods, 1))])
G = co.sparse([ \
co.spmatrix(-1.0, range(NumberOfPeriods), \
range(NumberOfAssets, NumberOfElements)), \
co.sparse([ \
[co.matrix(-Deviation)], \
[co.spmatrix(-1.0, range(NumberOfPeriods), range(NumberOfPeriods))]\
]), \
co.spmatrix(-1.0, range(NumberOfAssets), range(NumberOfAssets), \
(NumberOfAssets, NumberOfElements)) \
])
h = co.matrix(0.0, (2*NumberOfPeriods + NumberOfAssets, 1))
A = co.sparse([ \
[co.matrix(np.c_[MeanVector, np.ones(NumberOfAssets)]).trans()], \
[co.spmatrix([], [], [], (2, NumberOfPeriods))] \
])
b = co.matrix([TargetReturn, 1.0])
return co.solvers.lp(c, G, h, A, b)
# 最適期待ショートフォール・ポートフォリオの計算
def ES_Portfolio(ReturnData, TargetReturn, Alpha):
## ReturnData: 収益率データ
## TargetReturn: 目標収益率
## Alpha: 下方確率
## Output: 最適ポートフォリオ
co.solvers.options['show_progress'] = False
MeanVector = np.mean(ReturnData, axis=0)
NumberOfPeriods = ReturnData.shape[0]
NumberOfAssets = ReturnData.shape[1]
NumberOfElements = NumberOfAssets + NumberOfPeriods
c = co.matrix( \
[co.matrix(0.0, (NumberOfAssets, 1)), \
co.matrix(1.0/(Alpha*NumberOfPeriods), (NumberOfPeriods, 1)), -1.0])
G = co.sparse([ \
co.spmatrix(-1.0, range(NumberOfPeriods), \
range(NumberOfAssets, NumberOfElements), \
(NumberOfPeriods, NumberOfElements + 1)), \
co.sparse([ \
[co.matrix(-ReturnData)], \
[co.spmatrix(-1.0, range(NumberOfPeriods), range(NumberOfPeriods))], \
[co.matrix(1.0, (NumberOfPeriods, 1))] \
]), \
co.spmatrix(-1.0, range(NumberOfAssets), range(NumberOfAssets), \
(NumberOfAssets, NumberOfElements + 1)) \
])
h = co.matrix(0.0, (2*NumberOfPeriods + NumberOfAssets, 1))
A = co.sparse([ \
[co.matrix(np.c_[MeanVector, np.ones(NumberOfAssets)]).trans()], \
[co.spmatrix([], [], [], (2, NumberOfPeriods + 1))] \
])
b = co.matrix([TargetReturn, 1.0])
return co.solvers.lp(c, G, h, A, b)
# 株価データの読み込み
stockvalue = pd.read_csv('stocks.csv', index_col=0)
# 最小分散フロンティアの計算
R = (stockvalue.diff()/stockvalue.shift(1))[1:] * 100.0
T = R.shape[0]
mu = R.mean().values
Sigma = (R.cov().values)*((T-1.0)/T)
V_Target = np.linspace(mu.min(), mu.max(), num=50)
V_SD = np.array( \
[np.sqrt(2.0*MV_Portfolio(R.values, Target)['primal objective']) \
for Target in V_Target])
V_AD = np.array( \
[2.0*AD_Portfolio(R.values, Target)['primal objective'] \
for Target in V_Target])
V_ES = np.array( \
[ES_Portfolio(R.values, Target, 0.1)['primal objective'] \
for Target in V_Target])
# フロンティアのグラフの作成
()
fig1 = plt.figure(1, facecolor='w')
plt.subplot(1, 3, 1)
plt.plot(V_SD, V_Target, 'b-')
plt.plot(np.sqrt(np.diagonal(Sigma)), mu, 'rx')
plt.xlabel(u'standard deviation')
plt.ylabel(u'Expected rate of return')
plt.subplot(1, 3, 2)
plt.plot(V_AD, V_Target, 'b-')
plt.plot((R - R.mean()).abs().mean().values, mu, 'rx')
plt.xlabel(u'Absolute deviation')
plt.subplot(1, 3, 3)
plt.plot(V_ES, V_Target, 'b-')
plt.plot((-R[R <= R.quantile(0.1)]).mean().values, mu, 'rx')
plt.xlabel(u'Expected shortfall')
[此贴子已经被作者于2017-2-5 20:11编辑过]