这个算是比较基本的概念了。之前一直思考多股票交易的策略如何分配仓位（weight）。

之前策略的做法都是比较粗糙的等权重买入。

投资组合理论似乎给出了一个比较满意的答案。

1. 找到有效前沿。在既定的收益率下使组合的方差最小。

2. 找到sharpe最优的组合（收益-风险均衡点）

3. 找到风险最小的组合

效果图

有效前沿在资产配置领域里往后的应用，是结合个人的风险效用曲线来找到合适的投资组合交点。这一步此研究没有涉及。

正态性检验和蒙特卡洛完成投资组合优化¶

by 陈小米。

最近一直在思考怎样有效的配置资产组合。 很多时候根据条件选好股票池之后，通常简单粗暴的等分仓位给每只股票。 其实，这个过程中有很多可以优化的空间。

下面，给大家分享一下如何运用有效前沿进行资产组合优化。

PART ONE: 正态性检验¶

这部分是附赠福利。只对资产组合优化感兴趣的朋友可以直接跳到PART TWO。

1.导入模块¶

import pandas as pdimport numpy as npimport statsmodels.api as smimport scipy.stats as scsimport matplotlib.pyplot as plt

2.选取几只感兴趣的股票。¶

002697 红旗连锁, 600783 鲁信创投, 000413 东旭光电, 601588 北辰实业

stock = ['002697.XSHE','600783.XSHG','000413.XSHE','601588.XSHG']start_date = '2015-01-01'end_date = '2015-12-31'df = get_price(stock, start_date, end_date, 'daily',['close'])data = df['close']data.head()

3.比较一下机制股票的情况。规范起点为100.¶

(data/data.ix[0]*100).plot(figsize = (8,6))

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

4.计算收益率¶

用pandas计算收益率会比Numpy效率高一些，可以用shift方法

log_returns = np.log(data/data.shift(1))log_returns.head()

log_returns.hist(bins = 50, figsize = (9,6))

array([[<matplotlib.axes._subplots.AxesSubplot object at 0x7fd0420a5d90>,
        <matplotlib.axes._subplots.AxesSubplot object at 0x7fd04201bc50>],
       [<matplotlib.axes._subplots.AxesSubplot object at 0x7fd041f9ebd0>,
        <matplotlib.axes._subplots.AxesSubplot object at 0x7fd041f03650>]], dtype=object)

从以上结果看到4个数据集都和正态分布的要求相差太多。

5.输出每只股票的统计数据¶

#定义print_statistics函数，为了更加易于理解的方式#输出给定(历史或者模拟)数据集均值、偏斜度或者峰度等统计数字def print_statistics(array):sta = scs.describe(array)print '%14s %15s' %('statistic','value')print 30*'-'print '%14s %15d' %('size', sta[0])print '%14s %15.5f' %('min', sta[1][0])print '%14s %15.5f' %('max', sta[1][1])print '%14s %15.5f' %('mean', sta[2])print '%14s %15.5f' %('std', np.sqrt(sta[3]))print '%14s %15.5f' %('skew', sta[4])print '%14s %15.5f' %('kurtosis', sta[5])for st in stock:print '\nResults for stock %s' %stprint 30*'-'log_data = np.array(log_returns[st].dropna())print_statistics(log_data)

Results for stock 002697.XSHE

     statistic           value

          size             243
           min        -0.10697
           max         0.09596
          mean         0.00196
           std         0.04652
          skew        -0.13239
      kurtosis         0.33900

Results for stock 600783.XSHG

     statistic           value

          size             243
           min        -0.10563
           max         0.09554
          mean         0.00144
           std         0.05250
          skew        -0.30149
      kurtosis        -0.30247

Results for stock 000413.XSHE

     statistic           value

          size             243
           min        -0.10569
           max         0.09651
          mean         0.00073
           std         0.04084
          skew        -0.20718
      kurtosis         0.78082

Results for stock 601588.XSHG

     statistic           value

          size             243
           min        -0.10629
           max         0.09660
          mean         0.00020
           std         0.04250
          skew        -0.46498
      kurtosis         0.80946

6.画qq图观察数据¶

下面是002697.XSHE 对数收益率 分位数-分位数图

sm.qqplot(log_returns['002697.XSHE'].dropna(),line = 's')plt.grid(True)plt.xlabel('theoretical quantiles')plt.ylabel('sample quantiles')

<matplotlib.text.Text at 0x7fd041b9ff90>

很显然，样本的分位数值不在一条直线上，表明“非正态性”。左侧和右侧分别有许多值远低于和远高于直线。这是典型的Fat tails。 Fat tails是频数分布中观察到的两端的异常值。

7.进行正态性检验¶

def normality_test(array):'''    对给定的数据集进行正态性检验    组合了3中统计学测试    偏度测试（Skewtest）——足够接近0    峰度测试（Kurtosistest)——足够接近0    正态性测试    '''print 'Skew of data set %15.3f' % scs.skew(array)print 'Skew test p-value %14.3f' % scs.skewtest(array)[1]print 'Kurt of data set %15.3f' % scs.kurtosis(array)print 'Kurt test p-value %14.3f' % scs.kurtosistest(array)[1]print 'Norm test p-value %14.3f' % scs.normaltest(array)[1]for st in stock:print '\nResults for st %s' %stprint 32*'-'log_data = np.array(log_returns[st].dropna())normality_test(log_data)

Results for st 002697.XSHE

Skew of data set          -0.132
Skew test p-value          0.388
Kurt of data set           0.339
Kurt test p-value          0.233
Norm test p-value          0.338

Results for st 600783.XSHG

Skew of data set          -0.301
Skew test p-value          0.053
Kurt of data set          -0.302
Kurt test p-value          0.346
Norm test p-value          0.099

Results for st 000413.XSHE

Skew of data set          -0.207
Skew test p-value          0.179
Kurt of data set           0.781
Kurt test p-value          0.032
Norm test p-value          0.041

Results for st 601588.XSHG

Skew of data set          -0.465
Skew test p-value          0.004
Kurt of data set           0.809
Kurt test p-value          0.028
Norm test p-value          0.001

从上述测试的p值来看，否定了数据集呈正态分布的测试假设。 这说明，股票市场收益率的正态假设不成立。

PART TWO：均值-方差投资组合理论¶

该理论基于用均值和方差来表述组合的优劣的前提。将选取几只股票，用蒙特卡洛模拟初步探究组合的有效前沿。

通过最大Sharpe和最小方差两种优化来找到最优的资产组合配置权重参数。

最后，刻画出可能的分布，两种最优以及组合的有效前沿。

1.选取几只感兴趣的股票¶

000413 东旭光电，000063 中兴通讯，002007 华兰生物，000001 平安银行，000002 万科A

并比较一下数据（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, end_date, 'daily', ['close'])data = df['close']#规范化后时序数据(data/data.ix[0]*100).plot(figsize = (8,5))

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

2.计算不同证券的均值、协方差¶

每年252个交易日，用每日收益得到年化收益。

计算投资资产的协方差是构建资产组合过程的核心部分。运用pandas内置方法生产协方差矩阵。

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

000413.XSHE    0.184516
000063.XSHE    0.176790
002007.XSHE    0.309077
000001.XSHE   -0.102059
000002.XSHE    0.547441
dtype: float64

returns.cov()*252

3.给不同资产随机分配初始权重¶

由于A股不允许建立空头头寸，所有的权重系数均在0-1之间

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

array([ 0.37505798,  0.21652754,  0.31590981,  0.06087709,  0.03162758])

4.计算预期组合年化收益、组合方差和组合标准差¶

np.sum(returns.mean()*weights)*252

0.21622558669017816

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

0.23595133640121463

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

0.4857482232609962

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.04plt.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 instance at 0x7fd04155e638>

6.投资组合优化1——sharpe最大¶

建立statistics函数来记录重要的投资组合统计数据（收益，方差和夏普比）

通过对约束最优问题的求解，得到最优解。其中约束是权重总和为1。

def statistics(weights):weights = np.array(weights)port_returns = np.sum(returns.mean()*weights)*252port_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

  status: 0
 success: True
    njev: 4
    nfev: 28
     fun: -1.1623048291871221
       x: array([ -3.60840218e-16,   2.24626781e-16,   1.63619563e-01,
        -2.27085639e-16,   8.36380437e-01])
 message: 'Optimization terminated successfully.'
     jac: array([  1.81575805e-01,   5.40387481e-01,   8.18073750e-05,
         1.03137662e+00,  -1.60038471e-05,   0.00000000e+00])
     nit: 4

得到的最优组合权重向量为：

opts['x'].round(3)

array([-0.   ,  0.   ,  0.164, -0.   ,  0.836])

sharpe最大的组合3个统计数据分别为：

#预期收益率、预期波动率、最优夏普指数statistics(opts['x']).round(3)

array([ 0.508,  0.437,  1.162])

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

  status: 0
 success: True
    njev: 7
    nfev: 50
     fun: 0.38542969450547221
       x: array([  1.14787640e-01,   3.28089742e-17,   2.09584008e-01,
         3.53487044e-01,   3.22141307e-01])
 message: 'Optimization terminated successfully.'
     jac: array([ 0.3851725 ,  0.43591119,  0.3861807 ,  0.3849672 ,  0.38553924,  0.        ])
     nit: 7

方差最小的最优组合权重向量及组合的统计数据分别为：

optv['x'].round(3)

array([ 0.115,  0.   ,  0.21 ,  0.353,  0.322])

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

array([ 0.226,  0.385,  0.587])

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 instance at 0x7fd040d72518>