Tag Archives: finance

[Python] 알고리즘 트레이딩 (2) – Python Zipline 으로 일단 주식 사보기

앞선 포스팅에서는 주식을 하는 사람이라면 다 아는 이동평균선을 Python으로 만들면서 장난을 쳐보았다. 이번 포스팅에서는 조금 더 진지하게 알고리즘 트레이딩을 구현해보려고 한다.

zipline은 Python으로 알고리즘 트레이딩을 할 수 있는 라이브러리다.

zipline은 마찬가지로 Python으로 짜여진 quantopian 이라고 하는 실시간 트레이딩 엔진을 사용하고 있다.

zipline 뿐만 아니라 quantopian 사이트에 방문해보는 것도 재밌을 것 같다.

관심 있는 분들은 꼭 한번 들어가 보길 바란다.

 

주저리 주저리 설명을 보는 것 보다, 실제로 작동하는 코드를 설명과 함께 보는 것이 더 효율적일 것 같다.

 

import zipline as zp
import zipline.utils.factory as zpf
from zipline.api import order, sid, symbol

import pandas as pd
import pandas.io.data as web
import numpy as np

from datetime import datetime
import matplotlib.pyplot as plt

class BuyGoogle(zp.TradingAlgorithm):
    
    trace = False
    
    def __init__(self, trace=False):
        BuyGoogle.trace = trace
        super(BuyGoogle, self).__init__()
        
    def initialize(context):
        print("Initialize")
        print(context)
        print("finish initialize")
            
    def handle_data(self, context):
        print("Handle Data...")
        print(context)
        self.order(symbol("GOOGL"), 1)
        print("Finish handling")

 

import 한 클래스들. 그리고 zipline의 trading algorithms 을 변수로 받아서 Google을 사는 클래스를 구현해보았다.

보면 알겠지만, 일단 그냥 아무런 의미 없이 한주를 사는 클래스다. 그 외에는 그냥 로그를 찍는 느낌으로 context 를 찍어낸다.

handle_Data를 상속받아서 작성했는데, 저 handle_Data는 하루에 한번 실행된다.

context 내의 내용을 받아서 알고리즘 트레이딩을 구현할 수 있을 것이다. 일단 그냥 나는 한주만 산다.

data = zpf.load_from_yahoo(stocks=['GOOGL'], 
                           indexes={}, 
                           start=datetime(1990, 1, 1), 
                           end=datetime(2015, 12, 31),
                           adjusted=False)
data.plot(figsize=(12,8))

GOOGL_1990~2015

pandas와 비슷한 방식으로 google의 종가를 가져올 수 있다.

 

이제 주식을 사보자.

# 2015년 1월 1일 부터 15일간 사본다.
result = BuyGoogle().run(data['2015-01-01':'2015-01-15'])

결과

Initialize
BuyGoogle(
    capital_base=100000.0
    sim_params=
SimulationParameters(
    period_start=2015-01-02 00:00:00+00:00,
    period_end=2015-01-15 00:00:00+00:00,
    capital_base=100000.0,
    data_frequency=daily,
    emission_rate=daily,
    first_open=2015-01-02 14:31:00+00:00,
    last_close=2015-01-15 21:00:00+00:00),
    initialized=False,
    slippage=VolumeShareSlippage(
    volume_limit=0.25,
    price_impact=0.1),
    commission=PerShare(cost=0.03, min trade cost=None),
    blotter=Blotter(
    transact_partial=(VolumeShareSlippage(
    volume_limit=0.25,
    price_impact=0.1), PerShare(cost=0.03, min trade cost=None)),
    open_orders=defaultdict(<class 'list'>, {}),
    orders={},
    new_orders=[],
    current_dt=2006-01-01 00:00:00+00:00),
    recorded_vars={})
finish...
Handle Data...
BarData(SortedDict(<function AlgorithmSimulator.__init__.<locals>._get_removal_date at 0x7ff6904ebe18>, [(0, SIDData({'price': 529.549988, 'type': 4, '_sid': 0, '_freqstr': None, 'volume': 1000000000, 'sid': 0, '_initial_len': 3, 'dt': Timestamp('2015-01-02 00:00:00+0000', tz='UTC'), 'source_id': 'DataFrameSource-f0269ab25f28863c365bfd07c474cf5b'}))]))
Finish.
...........생략....................
Handle Data...
BarData(SortedDict(<function AlgorithmSimulator.__init__.<locals>._get_removal_date at 0x7ff6904ebe18>, [(0, SIDData({'price': 504.01000999999997, 'type': 4, '_sid': 0, '_freqstr': None, 'volume': 1000000000, 'sid': 0, '_initial_len': 3, 'dt': Timestamp('2015-01-15 00:00:00+0000', tz='UTC'), 'source_id': 'DataFrameSource-f0269ab25f28863c365bfd07c474cf5b'}))]))
Finish.
result
 algo_volatility  algorithm_period_return     alpha  \
2015-01-02 21:00:00         0.000000             0.000000e+00 -0.021200   
2015-01-05 21:00:00         0.000003            -3.000000e-07 -0.020398   
2015-01-06 21:00:00         0.001176            -1.288001e-04 -0.030992   
2015-01-07 21:00:00         0.000966            -1.589005e-04 -0.029256   
2015-01-08 21:00:00         0.001062            -1.064002e-04 -0.025784   
2015-01-09 21:00:00         0.001749            -3.543003e-04 -0.033546   
2015-01-12 21:00:00         0.001763            -5.376005e-04 -0.036167   
2015-01-13 21:00:00         0.002602            -2.535011e-04 -0.024737   
2015-01-14 21:00:00         0.002967             3.529929e-05 -0.015598   
2015-01-15 21:00:00         0.002907            -1.185993e-04 -0.017516   

                     benchmark_period_return           ...            \
2015-01-02 21:00:00                -0.000340           ...             
2015-01-05 21:00:00                -0.018612           ...             
2015-01-06 21:00:00                -0.027340           ...             
2015-01-07 21:00:00                -0.016028           ...             
2015-01-08 21:00:00                 0.001574           ...             
2015-01-09 21:00:00                -0.006843           ...             
2015-01-12 21:00:00                -0.014882           ...             
2015-01-13 21:00:00                -0.017422           ...             
2015-01-14 21:00:00                -0.023134           ...             
2015-01-15 21:00:00                -0.032168           ...             

                     starting_value  trading_days  \
2015-01-02 21:00:00        0.000000             1   
2015-01-05 21:00:00        0.000000             2   
2015-01-06 21:00:00      519.460022             3   
2015-01-07 21:00:00     1013.280030             4   
2015-01-08 21:00:00     1515.449982             5   
2015-01-09 21:00:00     2027.640016             6   
2015-01-12 21:00:00     2503.600005             7   
2015-01-13 21:00:00     2982.359988             8   
2015-01-14 21:00:00     3512.599916             9   
2015-01-15 21:00:00     4047.439944            10   

                                                          transactions  \
2015-01-02 21:00:00                                                 []   
2015-01-05 21:00:00  [{'amount': 1, 'price': 519.490022, 'order_id'...   
2015-01-06 21:00:00  [{'amount': 1, 'price': 506.67001500000003, 'o...   
2015-01-07 21:00:00  [{'amount': 1, 'price': 505.179994, 'order_id'...   
2015-01-08 21:00:00  [{'amount': 1, 'price': 506.940004, 'order_id'...   
2015-01-09 21:00:00  [{'amount': 1, 'price': 500.75000099999994, 'o...   
2015-01-12 21:00:00  [{'amount': 1, 'price': 497.08999800000004, 'o...   
2015-01-13 21:00:00  [{'amount': 1, 'price': 501.82998799999996, 'o...   
2015-01-14 21:00:00  [{'amount': 1, 'price': 505.95999299999994, 'o...   
2015-01-15 21:00:00  [{'amount': 1, 'price': 504.04000999999994, 'o...   

                     treasury_period_return  
2015-01-02 21:00:00                  0.0212  
2015-01-05 21:00:00                  0.0204  
2015-01-06 21:00:00                  0.0197  
2015-01-07 21:00:00                  0.0196  
2015-01-08 21:00:00                  0.0203  
2015-01-09 21:00:00                  0.0198  
2015-01-12 21:00:00                  0.0192  
2015-01-13 21:00:00                  0.0191  
2015-01-14 21:00:00                  0.0186  
2015-01-15 21:00:00                  0.0177  

[10 rows x 38 columns]

채결했던 주문에 대해서도 볼 수 있다.

result.iloc[1].orders
[{'amount': 1,
  'commission': 0.03,
  'created': Timestamp('2015-01-02 00:00:00+0000', tz='UTC'),
  'dt': Timestamp('2015-01-05 00:00:00+0000', tz='UTC'),
  'filled': 1,
  'id': 'b540e952e5b24d42ba533b49b1460bc0',
  'limit': None,
  'limit_reached': False,
  'reason': None,
  'sid': Equity(0, symbol='GOOGL', asset_name=None, exchange=None, start_date=Timestamp('1970-01-01 00:00:00+0000', tz='UTC'), end_date=Timestamp('2116-02-20 23:53:38.427387903+0000', tz='UTC'), first_traded=None, auto_close_date=None),
  'status': 1,
  'stop': None,
  'stop_reached': False},
 {'amount': 1,
  'commission': None,
  'created': Timestamp('2015-01-05 00:00:00+0000', tz='UTC'),
  'dt': Timestamp('2015-01-05 00:00:00+0000', tz='UTC'),
  'filled': 0,
  'id': 'af6837dcf91f41dfb8dca39675fea87f',
  'limit': None,
  'limit_reached': False,
  'reason': None,
  'sid': Equity(0, symbol='GOOGL', asset_name=None, exchange=None, start_date=Timestamp('1970-01-01 00:00:00+0000', tz='UTC'), end_date=Timestamp('2116-02-20 23:53:38.427387903+0000', tz='UTC'), first_traded=None, auto_close_date=None),
  'status': 0,
  'stop': None,
  'stop_reached': False}]

그래서 결과는?

result[['starting_cash', 'ending_cash', 'ending_value']]
 starting_cash    ending_cash  ending_value
2015-01-02 21:00:00  100000.000000  100000.000000      0.000000
2015-01-05 21:00:00  100000.000000   99480.509978    519.460022
2015-01-06 21:00:00   99480.509978   98973.839963   1013.280030
2015-01-07 21:00:00   98973.839963   98468.659969   1515.449982
2015-01-08 21:00:00   98468.659969   97961.719965   2027.640016
2015-01-09 21:00:00   97961.719965   97460.969964   2503.600005
2015-01-12 21:00:00   97460.969964   96963.879966   2982.359988
2015-01-13 21:00:00   96963.879966   96462.049978   3512.599916
2015-01-14 21:00:00   96462.049978   95956.089985   4047.439944
2015-01-15 21:00:00   95956.089985   95452.049975   4536.090090

왠지 약간 손해본 느낌이다.

좀 더 정확히 보자

result.portfolio_value
2015-01-02 21:00:00             NaN
2015-01-05 21:00:00   -3.000000e-07
2015-01-06 21:00:00   -1.285001e-04
2015-01-07 21:00:00   -3.010430e-05
2015-01-08 21:00:00    5.250864e-05
2015-01-09 21:00:00   -2.479265e-04
2015-01-12 21:00:00   -1.833651e-04
2015-01-13 21:00:00    2.842522e-04
2015-01-14 21:00:00    2.888736e-04
2015-01-15 21:00:00   -1.538932e-04
Name: portfolio_value, dtype: float64
result['returns']
2015-01-02 21:00:00    0.000000e+00
2015-01-05 21:00:00   -3.000000e-07
2015-01-06 21:00:00   -1.285001e-04
2015-01-07 21:00:00   -3.010430e-05
2015-01-08 21:00:00    5.250864e-05
2015-01-09 21:00:00   -2.479265e-04
2015-01-12 21:00:00   -1.833651e-04
2015-01-13 21:00:00    2.842522e-04
2015-01-14 21:00:00    2.888736e-04
2015-01-15 21:00:00   -1.538932e-04
Name: returns, dtype: float64

그렇다. 손해를 봤다.

result.portfolio_value.plot(figsize=(12, 8))

resut_GOOGL_201501

 

좀 더 화끈하게, 2015년 내내 한주씩 사보자.

result_for_2015 = BuyGoogle().run(data['2015'])
result_for_2015.portfolio_value.plot(figsize=(12, 8))

GOOGL_2015

외쳐 Google!

 

[Pandas] pandas.io.data.datareader를 이용한 주가 분석

pandas.io.data.datareader 는 정말 미친기능이다.

여기를 보면 알겠지만, 사람들이 자주 가져오는 주식등의 데이터를 순식간에 가져온다.

이렇게 지원하는데, 이것만 해도 왠만한 걸 다 할 수 있다.

<Parameter>

name  – String: 데이터셋의 Name이다. Google Finance의 경우 Stock의 symbol (예를 들어 Nasdaq Composite은 .IXIC다) 을 넣으면 된다.

data_sorce – String: ‘yahoo’, ‘google’ 등. 자세한건 여기 참조

start – String: start date다. ‘2015-01-01’ 이런식으로 넣으면 된다. Default 값은 ‘2010-1-1’ 이다.

end – String: end date. Default 값은 오늘이다.

 

자 한번 활용을 해보자.

import pandas.io.data as web
%matplotlib inline
DAX = web.DataReader(name='^GDAXI', data_source='yahoo', start='2000-1-1') # Yahoo에서 german DAX Index를 가져온다.

DAX.tail()
Open High Low Close Volume Adj Close
Date
2015-12-22 10598.190430 10623.980469 10400.610352 10488.750000 59880200 10488.750000
2015-12-23 10623.559570 10743.320312 10595.269531 10727.639648 69265000 10727.639648
2015-12-28 10748.370117 10756.169922 10627.459961 10653.910156 34275800 10653.910156
2015-12-29 10744.959961 10860.139648 10731.629883 10860.139648 51747000 10860.139648
2015-12-30 10855.169922 10857.429688 10743.009766 10743.009766 34183100 10743.009766

 

DAX['Close'].plot(figsize=(8, 5))

dax1

몇줄의 코드 만으로 순식간에 주식정보를 뽑아올 뿐만 아니라, 그래프까지 정확하게 그릴 수 있다.

 

자 그럼, 이번에는 빈 컬럼을 하나 생성해서 전날 주가와 비교한 로그값을 꽂아 보도록 하자.

%%time

import numpy as np

DAX['Ret_Loop'] = 0.0
for i in range(1, len(DAX)):
    DAX['Ret_Loop'][i] = np.log(DAX['Close'][i] / DAX['Close'][i - 1])

나의 t2.micro가 또 뻗은줄 알았다. 진지하게 성능 업을 고려해봐야겠다.

CPU times: user 1min 41s, sys: 0 ns, total: 1min 41s
Wall time: 1min 42s

DAX[['Close', 'Ret_Loop']].tail()
Close Ret_Loop
Date
2015-12-22 10488.750000 -0.000860
2015-12-23 10727.639648 0.022520
2015-12-28 10653.910156 -0.006897
2015-12-29 10860.139648 0.019172
2015-12-30 10743.009766 -0.010844

 

또는 이렇게 미친놈처럼  For문을 돌지 않고 다음과 같은 방법으로도 할 수 있다.

%time DAX['Return'] = np.log(DAX['Close'] / DAX['Close'].shift(1))

CPU times: user 1.82 ms, sys: 0 ns, total: 1.82 ms
Wall time: 3.01 ms

For문을 돌지 않았기 때문에 속도도 빠르고 깔끔하다. 결과는 같기 때문에 의심하지말자.

# 같은 값은 지우자.
del DAX['Ret_Loop']
DAX[['Close', 'Return']].plot(subplots=True, style='b', figsize=(20, 5))

dax2

주가 변동성이 옵션트레이더들에게 중요한 요소라고 한다면,
주식을 하는 사람들에게는 moving average, 즉 trend가 중요할 것이다.
moving average는 요즘쓰이는 기술지표들 중에 가장 유명하고 자주 쓰이는 요소다. 링크에 아주 잘 설명이 되어 있다.
이전 포스팅에서도 언급했지만, 아주 쉽게 만들 수 있다.

import pandas as pd
DAX['42d'] = pd.rolling_mean(DAX['Close'], window=42)
DAX['252d'] = pd.rolling_mean(DAX['Close'], window=252)
DAX[['Close', '42d', '252d']].plot(figsize=(20, 10))

dax3
옵션 트레이더에게는, 아래와 같은 moving historical volatility 정보가 유용할 지도 모르겠다.

import math
DAX['Mov_vol'] = pd.rolling_std(DAX['Return'], window=252) * math.sqrt(252)
DAX[['Close', 'Mov_vol', 'Return']].plot(subplots=True, style='b', figsize=(20, 10))

dax4

[Python] 주식

Python의 다양한 도구들을 이용하면, 주가를 쉽게 가져오고 또 원하는 데이터를 순식간에 뽑아낼 수 있다.
이번엔 주식, 주가를 가지고 놀아보도록 하겠당.

import numpy as np
import pandas as pd
import pandas.io.data as web

# yahoo에서 sp500 지수를 가져온다.
sp500 = web.DataReader('^GSPC', data_source='yahoo', start='1/1/2000', end='31/12/2015')

%matplotlib inline
sp500['Close'].plot(grid=True, figsize=(20, 10))

stock1

단기트렌드(42일)과 장기트렌드(252일)을 비교 해보자.

#42-day Trend
sp500['42d'] = np.round(pd.rolling_mean(sp500['Close'], window=42), 2)
#252-day Trend
sp500['252d'] = np.round(pd.rolling_mean(sp500['Close'], window=252), 2)

sp500[['Close', '42d', '252d']].tail()

 

Close 42d 252d
Date
2015-12-23 2064.290039 2068.93 2061.45
2015-12-24 2060.989990 2068.69 2061.37
2015-12-28 2056.500000 2068.47 2061.24
2015-12-29 2078.360107 2068.18 2061.19
2015-12-30 2063.360107 2067.56 2061.13
sp500[['Close', '42d', '252d']].plot(grid=True, figsize=(20, 10))

stock2

나름대로의 Trending Signal을 분석해본다고 가정하자.

Buy Signal: 42일 지표가 252일 지표 위에 있다.
중립: 42일 지표와 252일 지표에 큰차이가 없다.
Sell Signal: 42일 지표가 252일 지표 아래에 있다.

위와 같은 짓을 하기위해서는, 두 지표의 차이를 계산해야 할 것이다.

sp500['42-252'] = sp500['42d'] - sp500['252d']

그리고 50을 한계점이라고 가정하고, Signal을 판단해보자.+

SD = 50
sp500['Regime'] = np.where(sp500['42-252'] > SD, 1, 0)
sp500['Regime'] = np.where(sp500['42-252'] < -SD, -1, sp500['Regime'])
sp500['Regime'].plot(lw=1.5)
plt.ylim([-1.1, 1.1])

stock3

sp500['Market'] = np.log(sp500['Close'] / sp500['Close'].shift(1))
sp500['Strategy'] = sp500['Regime'].shift(1) * sp500['Market']
sp500[['Market', 'Strategy']].cumsum().apply(np.exp).plot(grid=True, figsize=(8, 5))

stock4

시장가격을 보면, 일단 경기침체가 있었던 2003, 2008~9기간에 short 전략을 가져갔던 사람들이 이익을 봤을 것으로 보인다.

그리고 이 기간을 제외하면, 가정했던 전략이 꽤 시장상황과 비슷한 추세를 보이는 것을 알 수 있다.

물론, 실제 투자에 필요한 몇몇 변수들이 완전히 무시되어서 그다지 쓸모있는 정보는 아닐 것 같다.

[Python] Monte Carlo Estimator

여기 와 연동되서 계속 하는 거다.

일단 소스 옮기는게 급해서 주석없이 복사한다.

 

from time import time
from math import exp, sqrt, log
from random import gauss, seed

seed(20000)
t0 = time()

S0 = 100. # Initial value
K  = 105. # Strike value
T  = 1.0  # maturity
r  = 0.05 # riskless short rate
sigma = 0.2 # volatility
M  = 50 # number of time stpes
dt = T / M #length of time interveal
I = 250000 # number of paths

S = []
for i in range(I):
    path = []
    for t in range(M + 1):
        if t == 0:
            path.append(S0)
        else:
            z = gauss(0.0, 1.0)
            St = path[t - 1] * exp((r - 0.5 * sigma ** 2) * dt + sigma * sqrt(dt) * z)
            path.append(St)
    S.append(path)

# Monte Carlo Estimator
C0 = exp(-r * T) * sum([max(path[-1] - K, 0) for path in S]) / I

tpy = time() - t0

print("European Option Value %7.3f'" % C0)
print("Dutaion in secnonds %7.3f'" % tpy)

log 버젼

import math
from numpy import *
from time import time

random.seed(20000)
t0 = time()

S0 = 100.
K = 105.
T = 1.0
r = 0.05
sigma = 0.2
M = 50
dt = T / M
I = 250000

S = S0 * exp(cumsum((r - 0.5 * sigma ** 2) * dt + sigma * math.sqrt(dt) * random.standard_normal((M+1, I)), axis = 0))
S[0] = S0

C0 = math.exp(-r * T) * sum(maximum(S[-1] - K, 0)) / I

tnp2 = time() - t0

print("European Option Value %7.3f" % C0)
print("Duration in Seconds %7.3f" % tnp2)

import matplotlib.pyplot as plt
plt.plot(S[:, :10])
plt.grid(True)
plt.xlabel('time step')
plt.ylabel('index level')

plt.show()

monte_Carlo_call_option_log

plt.hist(S[-1], bins=50)
plt.grid(True)
plt.xlabel('index level')
plt.ylabel('frequency')

plt.show()

Histogram of all simulated end-of-period index level values

plt.hist(maximum(S[-1] - K, 0), bins=50)
plt.grid(True)
plt.xlabel('option inner value')
plt.ylabel('frequency')
plt.ylim(0, 50000)

plt.show()

Histogram of all simulated end-of-period option inner values

[Python] Black-scholes-merton option pricing formula

Black-Scholes formula를 참고하자.

해설은 시간나면 쓰자.

Black-sholes Formula.py

\begin{align} C(S, t) &= N(d_1)S - N(d_2) Ke^{-r(T - t)} \\ d_1 &= \frac{1}{\sigma\sqrt{T - t}}\left[\ln\left(\frac{S}{K}\right) + \left(r + \frac{\sigma^2}{2}\right)(T - t)\right] \\ d_2 &= d_1 - \sigma\sqrt{T - t} \\ \end{align}


# bsm_콜 옵션가 구하는 공식
def bsm_call_value(S0, K, T, r, sigma):
   '''
    Parameters
    ==========
    S0 : float Initial Stock/index level
    K  : float 행사가격
    T  : float 만기일 (연 단위)
    r  : float risk-free short rate
    sigma : float 주식수익률의 변동성
    
    Returns
    =======
    value : float
        유로 콜옵션의 현재가격
    '''

    from math import log, sqrt, exp
    from scipy import stats

    S0 = float(S0)
    d1 = (log(S0 / K) + (r + 0.5 * sigma ** 2) * T) / (sigma * sqrt(T))
    d2 = (log(S0 / K) + (r - 0.5 * sigma ** 2) * T) / (sigma * sqrt(T))
    value = (S0 * stats.norm.cdf(d1, 0.0, 1.0)
            - K * exp(-r * T) * stats.norm.cdf(d2, 0.0, 1.0))
      # stats.norm.cdf --&amp;amp;amp;gt; cumulative distribution function
      #                    for normal distribution
    return value

# Vega function


def bsm_vega(S0, K, T, r, sigma):
    ''' Vega of European option in BSM model.
    Vega: measures sensitivity to volatility.
    이건 잘 모르겠다. Vega 구하는 공식이라고 하니 그냥 그러려니 한다.
    Returns
    =======
    vega : float
        

    '''
    from math import log, sqrt
    from scipy import stats

    S0 = float(S0)
    d1 = (log(S0 / K) + (r + 0.5 * sigma ** 2) * T) / (sigma * sqrt(T))
    vega = S0 * stats.norm.pdf(d1, 0.0, 1.0) * sqrt(T)
    return vega

# Implied volatility function


def bsm_call_imp_vol(S0, K, T, r, C0, sigma_est, it=100):
    ''' Implied volatility of European call option in BSM model.
    
    Returns
    =======
    simga_est : float
        numerically estimated implied volatility
    '''
    for i in range(it):
        sigma_est -= ((bsm_call_value(S0, K, T, r, sigma_est) - C0) / bsm_vega(S0, K, T, r, sigma_est))
    return sigma_est

시작!

# 2014-3-31 Closing Value 산정
V0 = 17.6339
# risk-free short rate
r = 0.01
h5 = pd.HDFStore('ipython/source/vstoxx_data_31032014.h5', 'r')
future_data = h5['futures_data']
options_data = h5['options_data']
h5.close()
DATE EXP_YEAR EXP_MONTH PRICE MATURITY TTM
496 2014-03-31 2014 4 17.85 2014-04-18 0.049
497 2014-03-31 2014 5 19.55 2014-05-16 0.126
498 2014-03-31 2014 6 19.95 2014-06-20 0.222
499 2014-03-31 2014 7 20.40 2014-07-18 0.299
500 2014-03-31 2014 8 20.70 2014-08-15 0.375
501 2014-03-31 2014 9 20.95 2014-09-19 0.471
502 2014-03-31 2014 10 21.05 2014-10-17 0.548
503 2014-03-31 2014 11 21.25 2014-11-21 0.644
options_data.info()

<class ‘pandas.core.frame.DataFrame’>
Int64Index: 395 entries, 46170 to 46564
Data columns (total 9 columns):
DATE 395 non-null int64
EXP_YEAR 395 non-null int64
EXP_MONTH 395 non-null int64
TYPE 395 non-null object
STRIKE 395 non-null float64
PRICE 395 non-null float64
MATURITY 395 non-null int64
TTM 395 non-null float64
IMP_VOL 395 non-null float64
dtypes: float64(4), int64(4), object(1)
memory usage: 30.9+ KB

options_data[['DATE', 'MATURITY', 'TTM', 'STRIKE', 'PRICE']].head()
DATE MATURITY TTM STRIKE PRICE
46170 1396224000000000000 1397779200000000000 0.049 1 16.85
46171 1396224000000000000 1397779200000000000 0.049 2 15.85
46172 1396224000000000000 1397779200000000000 0.049 3 14.85
46173 1396224000000000000 1397779200000000000 0.049 4 13.85
46174 1396224000000000000 1397779200000000000 0.049 5 12.85

날짜가 integer 64로 이상하게 나온다. 이럴때는 pd.to_Datetime() 을 써주면 간단하게 해결된다.

tol = 0.5  # tolerance level for moneyness
for option in options_data.index:
    # iterating over all option quotes
    forward = future_data[future_data['MATURITY'] == \ options_data.loc[option]['MATURITY']]['PRICE'].values[0]
      # picking the right futures value
    if (forward * (1 - tol) < options_data.loc[option]['STRIKE']
                             < forward * (1 + tol)):
        # only for options with moneyness within tolerance
        imp_vol = bsm_call_imp_vol(
                V0,  # VSTOXX value 
                options_data.loc[option]['STRIKE'],
                options_data.loc[option]['TTM'],
                r,   # short rate
                options_data.loc[option]['PRICE'],
                sigma_est=2.,  # estimate for implied volatility
                it=100)
        options_data['IMP_VOL'].loc[option] = imp_vol


plot_data = options_data[options_data['IMP_VOL'] > 0]
maturities = sorted(set(options_data['MATURITY']))

이제 그래프 그리면 끝!

import matplotlib.pyplot as plt
%matplotlib inline
plt.figure(figsize=(8, 6))
for maturity in maturities:
    data = plot_data[options_data.MATURITY == maturity]
    
    plt.plot(data['STRIKE'], data['IMP_VOL'], label=maturity.date(), lw=1.5)
    plt.plot(data['STRIKE'], data['IMP_VOL'], 'r,')
plt.grid(True)
plt.xlabel('strike')
plt.ylabel('implied volatility of volatility')
plt.legend()

bms