```python
import numpy as np
import scipy.stats as si
import scipy
import sympy as sy
import matplotlib.pyplot as plt
import pandas as pd
# import sympy.statistics as systats
```


```python
def euro_opt(S, K, T, r, sigma, option = 'call'):
    
    #S: spot price
    #K: strike price
    #T: time to maturity
    #r: interest rate
    #sigma: volatility of underlying asset
    
    d1 = (np.log(S / K) + (r + 0.5 * sigma ** 2) * T) / (sigma * np.sqrt(T))
    d2 = (np.log(S / K) + (r - 0.5 * sigma ** 2) * T) / (sigma * np.sqrt(T))
    
    if option == 'call':
        result = (S * si.norm.cdf(d1, 0.0, 1.0) - K * np.exp(-r * T) * si.norm.cdf(d2, 0.0, 1.0))
    if option == 'put':
        result = (K * np.exp(-r * T) * si.norm.cdf(-d2, 0.0, 1.0) - S * si.norm.cdf(-d1, 0.0, 1.0))
        
    return result
```


```python
euro_opt(27,27.25,6,0.065,0.0000000000001)
```




    8.550200169925013




```python
def vega(S, K, T, r, sigma):

    #S: spot price
    #K: strike price
    #T: time to maturity
    #r: interest rate
    #sigma: volatility of underlying asset
    
    d1 = (np.log(S / K) + (r + 0.5 * sigma ** 2) * T) / (sigma * np.sqrt(T))
    
    vega = S * si.norm.cdf(d1, 0.0, 1.0) * np.sqrt(T)
    
    return vega
```


```python
vega(27,27.25,6,0.065,.001)
```




    66.1362230551458




```python
black_scholes(-1,27,48.5,0.016,5,.65,0)
```




    23.176669536768507




```python
#-------------------------------------------------------------------------------
# Name:        Scenario Analysis Black Scholes
# Purpose:
#
# Author:      Jamie
#
# Created:     20/06/2012
# Copyright:   (c) Jamie 2012
# Licence:     <your licence>
#-------------------------------------------------------------------------------

from scipy import stats
import math

def black_scholes (cp, s, k, t, rf, v):
        """ Price an option using the Black-Scholes model.
        s: initial stock price
        k: strike price
        t: expiration time
        v: volatility
        rf: risk-free rate
        cp: +1/-1 for call/put
        """

        d1 = (math.log(s/k)+(rf+0.5*math.pow(v,2))*t)/(v*math.sqrt(t))
        d2 = d1 - v*math.sqrt(t)

        optprice = (cp*s*math.exp(-div*t)*stats.norm.cdf(cp*d1)) - (cp*k*math.exp(-rf*t)*stats.norm.cdf(cp*d2))
        return optprice

```


```python

```




    23.176669536768507




```python
def bissecao(f,a,b, p, tol=1e-8):
    # YOUR CODE HERE
    # p = [cp, s, k, t, v, rf, div, price]
    count = 0
    while b-a >= tol:
        count += 1
        z = (a+b)/2
        if f(p[0],p[1],p[2],p[3],z,p[5],p[6],p[7]) <= tol:
            return z, count
        if f(p[0],p[1],p[2],p[3],a,p[5],p[6],p[7])*f(p[0],p[1],p[2],p[3],z,p[5],p[6],p[7]) < 0:
            b = z
        else:
            a = z
    return z, count
```


```python
def h(cp, s, k, t, v, rf, div, price):
    #  print(black_scholes(cp, s, k, t, v, rf, div) - price)
    return abs(black_scholes(cp, s, k, t, v, rf, div) - price)
```


```python
bissecao(h,0.1,5,[-1,27,48.5,0.016,0.1,0.065,0,21.34])
```

    0.26928565347123623
    0.1095862197097155
    0.26928565347123623
    1.0132549214567277
    0.26928565347123623
    1.0132549214567277
    1.581320063767489
    1.0132549214567277
    1.581320063767489
    1.9006520304136956
    1.581320063767489
    1.9006520304136956
    2.0676267546886606
    1.9006520304136956
    2.0676267546886606
    2.1527723919917783
    2.0676267546886606
    2.1527723919917783
    2.195739933640265
    2.1527723919917783
    2.195739933640265
    2.2173199815349633
    2.195739933640265
    2.2173199815349633
    2.2281337790224143
    2.2173199815349633
    2.2281337790224143
    2.233546584478237
    2.2281337790224143
    2.233546584478237
    2.236254459318289
    2.233546584478237
    2.236254459318289
    2.2376087641967324
    2.236254459318289
    2.2376087641967324
    2.2382860084294514
    2.2376087641967324
    2.2382860084294514
    2.2386246534852816
    2.2382860084294514
    2.2386246534852816
    2.2387939817469658
    2.2386246534852816
    2.2387939817469658
    2.2388786473111004
    2.2387939817469658
    2.2388786473111004
    2.2389209804514856
    2.2388786473111004
    2.2389209804514856
    2.2389421471112456
    2.2389209804514856
    2.2389421471112456
    2.2389527304635273
    2.2389421471112456
    2.2389527304635273
    2.238958022145262
    2.2389527304635273
    2.238958022145262
    2.238960667987527
    2.238958022145262
    2.238960667987527
    2.238961990909008
    2.238960667987527
    2.238961990909008
    2.238962652369832
    2.238961990909008
    2.238962652369832
    2.2389629831002793
    2.238962652369832
    2.2389629831002793
    2.2389631484654977
    2.2389629831002793
    2.2389631484654977
    2.2389632311481087
    2.2389631484654977
    2.2389632311481087
    2.238963272489425
    2.2389632311481087
    2.238963272489425
    2.238963293160076
    2.238963272489425
    2.238963293160076
    2.238963303495403
    2.238963293160076
    2.238963303495403





    (4.999999990873039, 29)




```python
black_scholes(-1,27,32,0.0164,2.55,0.065,0)
```




    6.801033221436356




```python
ts = np.linspace(0.1,5,200)
plt.plot(ts,[h(-1,27,32,0.0164,0.7,0.065,0,21.34)]*200)
```


```python
black_scholes(-1,27,32,0.0164,0.065,0.715)
```




    5.000239133333103




```python

```




    3900




```python
cp = -1
s = 27
k = 32
t = 0.016
rf = 0.065
price = 5
def h(vol):
#   print(black_scholes(cp, s, k, t, vol, rf, div) - price)
    return black_scholes(cp, s, k, t, rf, vol) - price

scipy.optimize.bisect(h,1e-6,5,xtol=1e-16)
```




    0.7199792850227894




```python

```


```python
cp = 1
s = 27
k = 26.25
t = 0.016
rf = 0.065
price = 1.04
def h(vol):
    return black_scholes(cp, s, k, t, rf, vol) - price
scipy.optimize.bisect(h,1e-6,5,xtol=1e-16)
```




    0.4237032042981779




```python
black_scholes(1, 27, 26.25, 0.016, 0.065, 0.42370)
```




    1.0399962926697413




```python
from scipy import stats
import math

def black_scholes (cp, s, k, t, rf, v):
        """ Price an option using the Black-Scholes model.
        s: initial stock price
        k: strike price
        t: expiration time
        v: volatility
        rf: risk-free rate
        cp: +1/-1 for call/put
        """

        d1 = (math.log(s/k)+(rf+0.5*math.pow(v,2))*t)/(v*math.sqrt(t))
        d2 = d1 - v*math.sqrt(t)

        optprice = (cp*s*math.exp(-div*t)*stats.norm.cdf(cp*d1)) - (cp*k*math.exp(-rf*t)*stats.norm.cdf(cp*d2))
        return optprice

def volat_impl(cp, s, k, t, rf, price):
    def h(vol):
        return black_scholes(cp, s, k, t, rf, vol) - price
    return scipy.optimize.bisect(h,1e-6,5,xtol=1e-16)
```


```python

```


```python
black_scholes(1,27,28,0.016,0.065,0.3914)
```




    0.18993518619015592




```python
data = pd.read_csv('BRFOODS.csv',';')
# data
```


```python
'''
Setting CONSTANTS
'''

sigla_acao = 'BRFS3'
empresa = 'BRFoods S.A.'
preco_acao = 27.00
dias = [6,28]


puts = data['Tipo'] == 'PUT'
calls = data['Tipo'] == 'CALL'
dia1706 = data['TF'] == '17-06-2019'
dia1507 = data['TF'] == '15-07-2019'
```


```python
## PUT Com 6 dias a vencer

## Buscando as informações no DataFrame
df_k = data[puts & dia1706].iloc[0:,2:3]
df_s = data[puts & dia1706].iloc[0:,3:4]
ks = df_k.values.flatten()
Ss = df_s.values.flatten()

## Setando o array com as volatilidades a serem plotadas
vs = []
for (k,s) in zip(ks,Ss):
    vs.append(volat_impl(-1,preco_acao,k,dias[0]/365,0.065,s))

## Plot do gráfico
plt.figure(figsize=(10,6))
plt.plot(ks,vs, marker='o')
plt.xlabel('Strikes')
plt.ylabel('Volatilidade')
plt.title('Gráfico Smile - Volatilidade x Strike - PUT Com 6 dias a vencer')
plt.grid()
plt.show()
```


```python
## CALL Com 6 dias a vencer

## Buscando as informações no DataFrame
df_k = data[calls & dia1706].iloc[0:,2:3]
df_s = data[calls & dia1706].iloc[0:,3:4]
ks = df_k.values.flatten()
Ss = df_s.values.flatten()

## Setando o array com as volatilidades a serem plotadas
vs = []
for (k,s) in zip(ks,Ss):
    vs.append(volat_impl(1,preco_acao,k,dias[0]/365,0.065,s))

## Plot do gráfico
plt.figure(figsize=(10,6))
plt.plot(ks,vs, marker='o')
plt.xlabel('Strikes')
plt.ylabel('Volatilidade')
plt.title('Gráfico Smile - Volatilidade x Strike - CALL Com 6 dias a vencer')
plt.grid()
plt.show()
```


```python
## PUT Com 28 dias a vencer

## Buscando as informações no DataFrame
df_k = data[puts & dia1507].iloc[0:,2:3]
df_s = data[puts & dia1507].iloc[0:,3:4]
ks = df_k.values.flatten()
Ss = df_s.values.flatten()

## Setando o array com as volatilidades a serem plotadas
vs = []
for (k,s) in zip(ks,Ss):
    vs.append(volat_impl(-1,preco_acao,k,dias[1]/365,0.065,s))

## Plot do gráfico
plt.figure(figsize=(10,6))
plt.plot(ks,vs, marker='o')
plt.xlabel('Strikes')
plt.ylabel('Volatilidade')
plt.title('Gráfico Smile - Volatilidade x Strike - PUT com 28 dias a vencer')
plt.grid()
plt.show()
```


```python
## CALL Com 28 dias a vencer

## Buscando as informações no DataFrame
df_k = data[calls & dia1507].iloc[0:,2:3]
df_s = data[calls & dia1507].iloc[0:,3:4]
ks = df_k.values.flatten()
Ss = df_s.values.flatten()

## Setando o array com as volatilidades a serem plotadas
vs = []
for (k,s) in zip(ks,Ss):
    vs.append(volat_impl(1,preco_acao,k,dias[1]/365,0.065,s))

## Plot do gráfico
plt.figure(figsize=(10,6))
plt.plot(ks,vs, marker='o')
plt.xlabel('Strikes')
plt.ylabel('Volatilidade')
plt.title('Gráfico Smile - Volatilidade x Strike - CALL com 28 dias a vencer')
plt.grid()
plt.show()
```

68e0d0a6a421c7230d4777e9ee66fc9ca6f595ae
ipynb
Jupyter Notebook
mod-mat-financas-I-2019-1/Project_3/proj3.ipynb
mirandagil/university-courses
e70ce5262555e84cffb13e53e139e7eec21e8907
MIT
2019-12-23T16:39:01.000Z
2019-12-23T16:39:01.000Z
mod-mat-financas-I-2019-1/Project_3/proj3.ipynb
mirandagil/university-courses
e70ce5262555e84cffb13e53e139e7eec21e8907
MIT
mod-mat-financas-I-2019-1/Project_3/proj3.ipynb
mirandagil/university-courses
e70ce5262555e84cffb13e53e139e7eec21e8907
MIT
Qwen/Qwen-72B
1. YES
2. YES