작성자: admin 작성일시: 2016-07-08 16:17:03 조회수: 463 다운로드: 34
카테고리: 금융 공학 태그목록:

포트폴리오 최적화

In [18]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

주가 데이터 수집

In [19]:
import pandas_datareader.data as wb

symbols = ['AAPL', 'MSFT', 'YHOO', 'DB', 'GLD']
noa = len(symbols)
data = pd.DataFrame()
for sym in symbols:
    data[sym] = wb.DataReader(sym, data_source='yahoo',
                              end='2014-09-12')['Adj Close']
data.columns = symbols
data.tail()
Out[19]:
AAPL MSFT YHOO DB GLD
Date
2014-09-08 95.195338 44.324409 41.810001 34.039067 120.730003
2014-09-09 94.837239 44.601017 40.779999 33.834365 120.870003
2014-09-10 97.750397 44.677325 41.139999 34.399735 120.260002
2014-09-11 98.166562 44.829937 41.259998 34.360744 119.470001
2014-09-12 98.389165 44.543790 42.880001 34.146292 118.379997
In [20]:
(data / data.ix[0] * 100).plot(figsize=(8, 5), grid=True)
plt.show()

수익률 계산

In [15]:
rets = np.log(data / data.shift(1))
rets.mean() * 252
Out[15]:
AAPL    0.267080
MSFT    0.114505
YHOO    0.196165
DB     -0.125174
GLD     0.016054
dtype: float64
In [16]:
rets.cov() * 252
Out[16]:
AAPL MSFT YHOO DB GLD
AAPL 0.072784 0.020459 0.023243 0.041027 0.005231
MSFT 0.020459 0.049402 0.024244 0.046089 0.002105
YHOO 0.023243 0.024244 0.093349 0.051538 -0.000864
DB 0.041027 0.046089 0.051538 0.177517 0.008777
GLD 0.005231 0.002105 -0.000864 0.008777 0.032406
In [22]:
sns.heatmap(rets.cov() * 252)
plt.show()

포트폴리오 이론

  • 각 자산의 수익률이 $r_i$이고 자산 비중이 $w_i$이면
  • 포트폴리오 수익률은 $$ \begin{eqnarray} \mu_p &=& \text{E} [r_p] \\ &=& \text{E} \left[ \sum_i w_i r_i \right] \\ &=& \sum_i w_i \text{E} [r_i] \\ &=& \sum_i w_i \mu_i \\ &=& w^T \mu \end{eqnarray} $$ 이 식에서 $\mu$ 는 자산의 기대 수익률 벡터
  • 포트폴리오 분산은 $$ \begin{eqnarray} \sigma_p^2 &=& \text{E} \left[ (r_p - \mu_p)^2 \right] \\ &=& \text{E} \left[ (w^T r - w^T \mu)^2 \right] \\ &=& \text{E} \left[ w^T(r-\mu)(r-\mu)^Tw \right] \\ &=& w^T \text{E} \left[ (r-\mu)(r-\mu)^T \right] w \\ &=& w^T \Sigma w \end{eqnarray} $$ 이 식에서 $\Sigma$ 는 자산 수익률의 공분산 행렬
In [28]:
np.random.seed(2)
weights = np.random.random(noa)
weights /= np.sum(weights)
weights
Out[28]:
array([ 0.23349275,  0.01388454,  0.2943663 ,  0.2331326 ,  0.22512381])
In [29]:
np.sum(rets.mean() * weights) * 252
Out[29]:
0.09612740192056152
In [27]:
np.dot(weights.T, np.dot(rets.cov() * 252, weights))
Out[27]:
0.040091365274925586
In [30]:
np.sqrt(np.dot(weights.T, np.dot(rets.cov() * 252, weights)))
Out[30]:
0.20022828290460262

포트폴리오 수익률 시뮬레이션

In [31]:
prets = []
pvols = []
for p in range (2500):
    weights = np.random.random(noa)
    weights /= np.sum(weights)
    prets.append(np.sum(rets.mean() * weights) * 252)
    pvols.append(np.sqrt(np.dot(weights.T, 
                        np.dot(rets.cov() * 252, weights))))
prets = np.array(prets)
pvols = np.array(pvols)
In [36]:
plt.scatter(pvols, prets, c=prets/pvols, marker='o', cmap=mpl.cm.jet)
plt.grid(True)
plt.xlabel('expected volatility')
plt.ylabel('expected return')
plt.colorbar(label='Sharpe ratio')
plt.show()

포트폴리오 통계

  • 수익률, 변동성, 샤프 지수 계산
In [37]:
def statistics(weights):
    ''' Return portfolio statistics.
    
    Parameters
    ==========
    weights : array-like
        포트폴리오 내의 증권 비중
    
    Returns
    =======
    pret : float
        포트폴리오 수익률의 기댓값
    pvol : float
        포트폴리오 변동성의 기댓값
    pret / pvol : float
        무위험 이자율이 0일 때의 포트폴리오 샤프 지수
    '''
    weights = np.array(weights)
    pret = np.sum(rets.mean() * weights) * 252
    pvol = np.sqrt(np.dot(weights.T, np.dot(rets.cov() * 252, weights)))
    return np.array([pret, pvol, pret / pvol])

샤프 지수 최대화

In [38]:
import scipy.optimize as sco
In [39]:
def min_func_sharpe(weights):
    return -statistics(weights)[2]
In [40]:
cons = ({'type': 'eq', 'fun': lambda x:  np.sum(x) - 1})
In [42]:
bnds = tuple((0, 1) for x in range(noa))
bnds
Out[42]:
((0, 1), (0, 1), (0, 1), (0, 1), (0, 1))
In [43]:
noa * [1. / noa,]
Out[43]:
[0.2, 0.2, 0.2, 0.2, 0.2]
In [44]:
%%time
opts = sco.minimize(min_func_sharpe, noa * [1. / noa,], method='SLSQP',
                       bounds=bnds, constraints=cons)
CPU times: user 40 ms, sys: 0 ns, total: 40 ms
Wall time: 39.6 ms
In [45]:
opts
Out[45]:
     fun: -1.0630084836641971
     jac: array([ -1.82956457e-04,  -7.02306628e-04,   7.18027353e-04,
         1.51409782e+00,   1.54867768e-03,   0.00000000e+00])
 message: 'Optimization terminated successfully.'
    nfev: 36
     nit: 5
    njev: 5
  status: 0
 success: True
       x: array([  6.61851751e-01,   8.64634623e-02,   2.51684786e-01,
         0.00000000e+00,   2.72744609e-17])
In [46]:
opts['x'].round(3)
Out[46]:
array([ 0.662,  0.086,  0.252,  0.   ,  0.   ])
In [47]:
statistics(opts['x']).round(3)
Out[47]:
array([ 0.236,  0.222,  1.063])
In [73]:
plt.scatter(pvols, prets, c=prets/pvols, marker='o', cmap=mpl.cm.jet)
plt.grid(True)
plt.xlabel('expected volatility')
plt.ylabel('expected return')
plt.colorbar(label='Sharpe ratio')

pt_opts = statistics(opts['x']).round(3)
plt.scatter(pt_opts[1], pt_opts[0], marker="*", s=500, alpha=0.5)
plt.show()