#引入需要的包
import pandas as pd
import numpy as np
import statsmodels.api as sm #统计运算
import scipy.stats as scs #科学计算
import matplotlib.pyplot as plt #绘图

#1.选取几只感兴趣的股票
#并比较一下数据（2015-01-01至2015-12-31）
stock_set = ['000413.XSHE','000063.XSHE','002007.XSHE','000001.XSHE','000002.XSHE']
noa = len(stock_set)
df = get_price(stock_set, start_date = '2015-01-01', end_date ='2015-12-31', frequency='daily', fields=['close'])
data = df['close']

print(str(data))
#规范化后时序数据
(data/data.ix[0]*100).plot(figsize = (8,5))

            000413.XSHE  000063.XSHE     ...       000001.XSHE  000002.XSHE
2015-01-05         7.40        15.44     ...             10.44        12.87
2015-01-06         7.40        16.34     ...             10.28        12.40
2015-01-07         7.40        16.01     ...             10.09        12.29
2015-01-08         7.40        16.24     ...              9.75        11.73
2015-01-09         7.40        15.88     ...              9.83        11.61
2015-01-12         7.40        15.55     ...              9.63        11.33
2015-01-13         7.40        15.61     ...              9.57        11.28
2015-01-14         7.40        15.57     ...              9.65        11.31
2015-01-15         7.40        15.96     ...             10.00        11.89
2015-01-16         7.40        16.06     ...             10.02        11.67
2015-01-19         7.40        14.62     ...              9.01        10.51
2015-01-20         7.40        15.50     ...              9.01        10.66
2015-01-21         7.40        15.99     ...              9.40        11.27
2015-01-22         7.40        16.01     ...              9.32        11.72
2015-01-23         7.40        15.99     ...              9.38        12.02
2015-01-26         7.40        16.56     ...              9.35        11.83
2015-01-27         7.40        16.32     ...              9.12        11.38
2015-01-28         8.14        15.99     ...              9.16        11.22
2015-01-29         8.03        15.70     ...              9.06        11.23
2015-01-30         7.73        15.58     ...              9.08        11.33
2015-02-02         7.58        15.16     ...              8.88        11.12
2015-02-03         7.67        15.54     ...              9.09        11.19
2015-02-04         7.57        15.14     ...              8.93        11.05
2015-02-05         7.35        15.23     ...              8.99        10.96
2015-02-06         7.20        14.82     ...              8.80        10.50
2015-02-09         7.41        14.87     ...              8.81        10.52
2015-02-10         7.39        15.23     ...              8.97        10.72
2015-02-11         7.43        16.00     ...              8.95        10.71
2015-02-12         7.45        15.91     ...              9.03        10.77
2015-02-13         7.56        16.00     ...              9.09        11.02
...                 ...          ...     ...               ...          ...
2015-11-20         9.84        18.57     ...              9.90        12.94
2015-11-23         9.48        18.13     ...              9.82        12.83
2015-11-24         9.73        18.70     ...              9.69        12.95
2015-11-25         9.72        18.72     ...              9.72        13.00
2015-11-26         9.83        18.48     ...              9.65        12.95
2015-11-27         9.12        17.28     ...              9.25        12.74
2015-11-30        10.04        17.18     ...              9.26        13.46
2015-12-01         9.97        17.39     ...              9.27        14.81
2015-12-02         9.65        17.87     ...              9.87        16.29
2015-12-03         9.77        17.83     ...              9.82        17.10
2015-12-04         9.58        17.43     ...              9.56        16.95
2015-12-07         9.73        17.69     ...              9.58        16.09
2015-12-08         9.25        17.20     ...              9.43        15.87
2015-12-09         8.93        17.22     ...              9.46        17.46
2015-12-10         8.86        17.36     ...              9.43        17.45
2015-12-11         8.90        17.38     ...              9.33        17.77
2015-12-14         9.06        17.79     ...              9.52        17.93
2015-12-15         9.07        18.38     ...              9.40        18.83
2015-12-16         9.12        18.04     ...              9.38        18.03
2015-12-17         9.35        18.32     ...              9.52        19.83
2015-12-18         9.17        18.46     ...              9.65        21.82
2015-12-21         9.36        19.08     ...              9.87        21.82
2015-12-22         9.39        18.79     ...              9.81        21.82
2015-12-23         9.14        18.48     ...              9.85        21.82
2015-12-24         9.13        18.58     ...              9.73        21.82
2015-12-25         9.09        18.69     ...              9.79        21.82
2015-12-28         8.79        18.05     ...              9.45        21.82
2015-12-29         8.92        18.30     ...              9.54        21.82
2015-12-30         9.04        18.83     ...              9.55        21.82
2015-12-31         8.83        18.32     ...              9.46        21.82

[244 rows x 5 columns]

/opt/conda/lib/python3.6/site-packages/ipykernel_launcher.py:10: DeprecationWarning: 
.ix is deprecated. Please use
.loc for label based indexing or
.iloc for positional indexing

See the documentation here:
http://pandas.pydata.org/pandas-docs/stable/indexing.html#ix-indexer-is-deprecated
  # Remove the CWD from sys.path while we load stuff.

<matplotlib.axes._subplots.AxesSubplot at 0x7f2f00b21a90>

#2.计算不同证券的均值、协方差
#每年252个交易日，用每日收益得到年化收益。计算投资资产的协方差是构建资产组合过程的核心部分。运用pandas内置方法生产协方差矩阵。

returns = np.log(data / data.shift(1))
returns.mean()*252

000413.XSHE    0.183219
000063.XSHE    0.177366
002007.XSHE    0.309165
000001.XSHE   -0.102223
000002.XSHE    0.547481
dtype: float64

returns.cov()*252

#3.给不同资产随机分配初始权重
#由于A股不允许建立空头头寸，所有的权重系数均在0-1之间

weights = np.random.random(noa)
weights /= np.sum(weights)

weights

array([0.20219897382356528, 0.24130281711887622, 0.0543513216091828,
       0.46183860222681006, 0.040308285221565654])

#4.计算预期组合年化收益、组合方差和组合标准差
np.sum(returns.mean()*weights)*252

0.07150660572458754

np.dot(weights.T, np.dot(returns.cov()*252,weights))

0.18008470701749085

np.sqrt(np.dot(weights.T, np.dot(returns.cov()* 252,weights)))

0.4243638851475121

'''
5.用蒙特卡洛模拟产生大量随机组合
进行到此，我们最想知道的是给定的一个股票池（证券组合）如何找到风险和收益平衡的位置。
下面通过一次蒙特卡洛模拟，产生大量随机的权重向量，并记录随机组合的预期收益和方差。
'''
port_returns = []
port_variance = []

for p in range(4000):

    weights = np.random.random(noa)
    weights /=np.sum(weights)

    port_returns.append(np.sum(returns.mean()*252*weights)) 
    port_variance.append(np.sqrt(np.dot(weights.T, np.dot(returns.cov()*252, weights))))

port_returns = np.array(port_returns)
port_variance = np.array(port_variance)
    
#无风险利率设定为4%
risk_free = 0.04

plt.figure(figsize = (8,4))
plt.scatter(port_variance, port_returns, c=(port_returns-risk_free)/port_variance, marker = 'o')
plt.grid(True)
plt.xlabel('excepted volatility')
plt.ylabel('expected return')
plt.colorbar(label = 'Sharpe ratio')

<matplotlib.colorbar.Colorbar at 0x7f2f00a03be0>

'''
6.投资组合优化1——sharpe最大
建立statistics函数来记录重要的投资组合统计数据（收益，方差和夏普比）
通过对约束最优问题的求解，得到最优解。其中约束是权重总和为1。
'''
def statistics(weights):
    weights = np.array(weights)
    port_returns = np.sum(returns.mean()*weights)*252
    port_variance = np.sqrt(np.dot(weights.T, np.dot(returns.cov()*252,weights)))

    return np.array([port_returns, port_variance, port_returns/port_variance])

#最优化投资组合的推导是一个约束最优化问题

import scipy.optimize as sco

#最小化夏普指数的负值

def min_sharpe(weights):
    return -statistics(weights)[2]

#约束是所有参数(权重)的总和为1。这可以用minimize函数的约定表达如下
cons = ({'type':'eq', 'fun':lambda x: np.sum(x)-1})

#我们还将参数值(权重)限制在0和1之间。这些值以多个元组组成的一个元组形式提供给最小化函数
bnds = tuple((0,1) for x in range(noa))

#优化函数调用中忽略的唯一输入是起始参数列表(对权重的初始猜测)。我们简单的使用平均分布。
opts = sco.minimize(min_sharpe, noa*[1./noa,], method = 'SLSQP', bounds = bnds, constraints = cons)

opts

     fun: -1.1634349336909844
     jac: array([0.1850493848323822, 0.5387051850557327, 8.337199687957764e-05,
       1.033293530344963, -1.633167266845703e-05])
 message: 'Optimization terminated successfully.'
    nfev: 28
     nit: 4
    njev: 4
  status: 0
 success: True
       x: array([0.0, 2.784665950614344e-16, 0.1638689231137425,
       3.6777449495635453e-16, 0.8361310768862571])

#得到的最优组合权重向量为：
opts['x'].round(3)

array([0.0, 0.0, 0.164, 0.0, 0.836])

#sharpe最大的组合3个统计数据分别为：
#预期收益率、预期波动率、最优夏普指数
statistics(opts['x']).round(3)

array([0.508, 0.437, 1.163])

'''
7.投资组合优化2——方差最小
接下来，我们通过方差最小来选出最优投资组合。
但是我们定义一个函数对 方差进行最小化
'''

def min_variance(weights):
    return statistics(weights)[1]

optv = sco.minimize(min_variance, noa*[1./noa,],method = 'SLSQP', bounds = bnds, constraints = cons)
optv

     fun: 0.38527596736539027
     jac: array([0.3850162401795387, 0.43552573025226593, 0.38600122556090355,
       0.3848104737699032, 0.38540058955550194])
 message: 'Optimization terminated successfully.'
    nfev: 50
     nit: 7
    njev: 7
  status: 0
 success: True
       x: array([0.11334405891423127, 3.581678817464809e-18, 0.2107763447259679,
       0.3523988464457304, 0.3234807499140705])

#方差最小的最优组合权重向量及组合的统计数据分别为：
optv['x'].round(3)

array([0.113, 0.0, 0.211, 0.352, 0.323])

#得到的预期收益率、波动率和夏普指数
statistics(optv['x']).round(3)

array([0.227, 0.385, 0.589])

'''
8.组合的有效前沿
有效前沿有既定的目标收益率下方差最小的投资组合构成。
在最优化时采用两个约束，1.给定目标收益率，2.投资组合权重和为1。
'''
def min_variance(weights):
    return statistics(weights)[1]

#在不同目标收益率水平（target_returns）循环时，最小化的一个约束条件会变化。
target_returns = np.linspace(0.0,0.5,50)
target_variance = []

for tar in target_returns:
    cons = ({'type':'eq','fun':lambda x:statistics(x)[0]-tar},{'type':'eq','fun':lambda x:np.sum(x)-1})
    res = sco.minimize(min_variance, noa*[1./noa,],method = 'SLSQP', bounds = bnds, constraints = cons)
    target_variance.append(res['fun'])

target_variance = np.array(target_variance)

'''
下面是最优化结果的展示。
叉号：构成的曲线是有效前沿（目标收益率下最优的投资组合）
红星：sharpe最大的投资组合
黄星：方差最小的投资组合
'''
plt.figure(figsize = (8,4))

#圆圈：蒙特卡洛随机产生的组合分布
plt.scatter(port_variance, port_returns, c = port_returns/port_variance,marker = 'o')

#叉号：有效前沿
plt.scatter(target_variance,target_returns, c = target_returns/target_variance, marker = 'x')

#红星：标记最高sharpe组合
plt.plot(statistics(opts['x'])[1], statistics(opts['x'])[0], 'r*', markersize = 15.0)

#黄星：标记最小方差组合
plt.plot(statistics(optv['x'])[1], statistics(optv['x'])[0], 'y*', markersize = 15.0)

plt.grid(True)
plt.xlabel('expected volatility')
plt.ylabel('expected return')
plt.colorbar(label = 'Sharpe ratio')

<matplotlib.colorbar.Colorbar at 0x7f2f009ad400>

	000413.XSHE	000063.XSHE	002007.XSHE	000001.XSHE	000002.XSHE
000413.XSHE	0.421851	0.206538	0.190642	0.096497	0.081273
000063.XSHE	0.206538	0.372291	0.212102	0.150713	0.143867
002007.XSHE	0.190642	0.212102	0.326678	0.096210	0.074879
000001.XSHE	0.096497	0.150713	0.096210	0.205699	0.138076
000002.XSHE	0.081273	0.143867	0.074879	0.138076	0.231268

量化交易吧 / 数理科学 帖子：3364470 新帖：36

用python实现Markowitz投资组合优化

外汇老老法师发表于：7 月 3 日 16：00回复(1)

全部回复

0/140

粉丝:733

帖子数:0

粉丝:555

帖子数:3

粉丝:555

帖子数:0

量化课程

热门标签

删除回复

确认要删除这篇文章么？

举报用户

信息提示

该文章已删除

设置置顶

完成设置【置顶】！

设置置顶

已取消设置【置顶】！

设置精华

完成设置【精华】！

设置精华

已取消设置【精华】！

审核信息

该文章已审核通过

审核信息

您已设置该文章审核不通过

举报成功

您已举报成功

用户登录

移动帖子

创建私信

屏蔽提示

确认要屏蔽该用户么？

屏蔽回复

您已对该用户实现屏蔽

信息回复

已发送成功

量化交易吧 / 数理科学帖子：3364470 新帖：36