# -*- coding: utf-8 -*- """ @author: Gabriel Turinici, 29/04/2020 """ import numpy as np import matplotlib.pyplot as plt from scipy.stats import norm def blsprice(Price,Strike,Rate,TimeToMaturity,Volatility,DividendRate=0): """input: S:Price - Current price of the underlying asset. Strike:Strike - Strike (i.e., exercise) price of the option. Rate: Rate - Annualized continuously compounded risk-free rate of return over the life of the option, expressed as a positive decimal number. TimeToMaturity:Time - Time to expiration of the option, expressed in years. Volatility: volatility DividendRate = continuous dividend rate """ if TimeToMaturity <= 1e-6: # the option already expired call = np.max(Price-Strike,0) put = np.max(Strike-Price,0) return call,put d1 = np.log(Price/Strike)+(Rate-DividendRate + Volatility**2/2.0)*TimeToMaturity; d1 = d1/(Volatility* np.sqrt(TimeToMaturity)) d2 = d1-(Volatility*np.sqrt(TimeToMaturity)) call = Price * np.exp(-DividendRate*TimeToMaturity) * norm.cdf(d1)-Strike* np.exp(-Rate*TimeToMaturity) * norm.cdf(d2) put = Strike* np.exp(-Rate*TimeToMaturity) * norm.cdf(-d2)-Price* np.exp(-DividendRate*TimeToMaturity) * norm.cdf(-d1) return call,put def blsdelta(Price,Strike,Rate,TimeToMaturity,Volatility,DividendRate=0): """input: S:Price - Current price of the underlying asset. Strike:Strike - Strike (i.e., exercise) price of the option. Rate: Rate - Annualized continuously compounded risk-free rate of return over the life of the option, expressed as a positive decimal number. TimeToMaturity:Time - Time to expiration of the option, expressed in years. Volatility: volatility DividendRate = continuous dividend rate """ if TimeToMaturity <= 1e-6: # the option already expired call = (Price>=Strike).astype(np.float) # 1 if in the money, zero otherwise put = (Price<=Strike).astype(np.float) # cf above return call,put d1 = np.log(Price/Strike)+(Rate-DividendRate + Volatility**2/2.0)*TimeToMaturity; d1 = d1/(Volatility* np.sqrt(TimeToMaturity)) call = np.exp(-DividendRate*TimeToMaturity) * norm.cdf(d1) put = -np.exp(-DividendRate*TimeToMaturity) * norm.cdf(-d1) return call,put T=1.0 N=255*8*4 M=1000 # nombre de scenarios dt = T/N W0=0 timerange = np.linspace(0,T,N+1,endpoint=True) timerange=timerange[:,None] dW= np.sqrt(dt)*np.random.randn(N,M) W=np.zeros((N+1,M)) W[0,:]=W0 W[1:,:]=W0+np.cumsum(dW,0) S0=100. mu=0.1 sigma=0.25 sigmareel=sigma*1.0 #sigmareel=sigma*0.90 #sigmareel=sigma*1.10 taux_r=0.05 def st(mu,sigma,t,wt): return np.exp((mu-sigma**2/2.)*t + sigma*wt) St = st(mu,sigmareel,timerange,W)*S0 plt.figure(5) for ii in range(15): plt.plot(timerange,St[:,ii]) K=110 prixcall,_ = blsprice(S0,K,taux_r,T,sigma) cash=np.zeros_like(St) deltaancien,_=blsdelta(S0,K,taux_r,T,sigma) #calcul de delta cash[0,:]=prixcall - deltaancien*S0 #prix du call qui a ete percu moins le #prix a payer pour la mise en place de la strategie for jj in range(N): cash[jj+1,:] =cash[jj,:]*np.exp(taux_r*dt)#actualisation deltanouveau,_=blsdelta(St[jj+1,:],K,taux_r, T-timerange[jj+1,0],sigma) #calcul du nouveau delta cash[jj+1,:]-= (deltanouveau- deltaancien)*St[jj+1,:] #paiement de l'actualisation de la #couverture en delta deltaancien=deltanouveau #livraison cash[-1,:] -= np.maximum(St[-1,:]-K,0) #livraison des gains cash[-1,:] += deltanouveau*St[-1,:]#vente portef de couverture plt.figure(6) plt.hist(cash[-1,:]*100/prixcall) plt.title('hist du cash - trading de vol- en pourcentages du prix ')