Dans la suite des exercices on supposera que l'actif risqué $S$ suit la dynamique suivante : $$ \begin{equation} S_t = S_0 exp \left( \left( r - \frac{\sigma^2}{2} \right) t + \sigma W_t \right) \end{equation} $$
Pour aller plus loin : Commençons par calculer une première fois le prix de l'option asiatique
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
import time
### paramètres "physique"
T_input = 1
r_input = 0.05
S0_input = 100
sigma_input = 0.2
K_input = 95
### paramètres "numériques"
n_input = 200*T_input
Mc_input = 10000
## Ici on construit la matrice de variance-covariance du mouvement brownien
def MatBrownien(T = T_input,n = n_input):
"""
Fonction créant la matrice de covariance du mouvement brownien. Soit {B_t}_{t \in [0;T]}, on alors \forall t,s \in [0;T],
Cov(B_t, B_s) = min(t,s). Comme de plus {B_t}_{t \in [0;T]} est un processus gaussien, il est caractérisé par
sa fonction de covariance.
--------
Arguments:
T : Temps final de la discrétisation (Maturité).
n : Nombre de points des discrétisations
Returns:
Mat : Matrice de covariance du mouvement brownien.
"""
t = np.linspace(0,T,n+1)
s = np.linspace(0,T,n+1)
T_mesh,S_mesh = np.meshgrid(t,s)
Mat = np.minimum(T_mesh,S_mesh)
return Mat
## Ici on crée la matrice W qui est censé représenter un mouvement brownien sur l'intervalle [0,T]
def dMouvementBrownien(T = T_input, n = n_input, d = Mc_input):
"""
Fonction permettant de créer la matrice W du mouvement brownien pour un échantillonage type Monte-Carlo.
On crée ici, d trajectoires browniennes, à partir d'une discrétisation de n+1 points de l'intervalle [0, T]
--------
Arguments:
T : Temps final de la discrétisation (Maturité).
n : Nombre de points de la discrétisation.
d : Nombre de trajectoires souhaitées. (Autant qu'il y a d'échantillonage de Monte-Carlo souhaité).
Returns:
W : Matrice des trajectoires browniennes. Taille : (d,n+1)
"""
t = np.linspace(0,T,n+1)
W = np.random.normal(size = (d,n+1))
Mat = MatBrownien(T,n)
A = np.linalg.cholesky(Mat[1:,1:])
B = np.zeros((n+1,n+1))
B[1:,1:] = A
for i in range(d):
W[i] = np.dot(B,W[i])
return (t,W)
Les fonctions ci-dessus nous permettent via la méthode de Cholesky, de générer d-mouvement browniens qui nous sont retournés sout la forme d'une matrice de MB : $$ \mathbb{G} = \begin{pmatrix} W_{t_{0}}^{1} & W_{t_{1}}^{1} & \cdots & W_{t_{n}}^{1}\\ W_{t_{0}}^{2} & W_{t_{1}}^{2} & \cdots & W_{t_{n}}^{2}\\ \vdots & \vdots & \ddots & \vdots \\ W_{t_{0}}^{d} & W_{t_{1}}^{d} & \cdots & W_{t_{n}}^{d}\\ \end{pmatrix} = \begin{pmatrix} 0 & W_{t_{1}}^{1} & \cdots & W_{t_{n}}^{1}\\ 0 & W_{t_{1}}^{2} & \cdots & W_{t_{n}}^{2}\\ \vdots & \vdots & \ddots & \vdots \\ 0 & W_{t_{1}}^{d} & \cdots & W_{t_{n}}^{d}\\ \end{pmatrix} $$
t,W = dMouvementBrownien(T_input,n_input,4)
plt.subplots(figsize=(8,5))
for i in range(4):
plt.plot(t,W[i], label = "Mvmt brownien n°{}".format(i+1))
plt.xlabel("temps")
plt.title("Mouvement brownien | n = {}".format(n_input))
plt.legend()
plt.show()
prixActif = lambda t,W,sigma,r, s0: s0*np.exp(sigma*W + t*(r - (sigma**2)/2 ))
En injectant la matrice de browniens calculé plus haut dans la fonction prixActif, on obtient la matrice suivante : $$ \text{Mat_{S_t}} = \begin{pmatrix} S_{t_{0}}^{1} & S_{t_{1}}^{1} & \cdots & S_{t_{n}}^{1}\\ S_{t_{0}}^{2} & S_{t_{1}}^{2} & \cdots & S_{t_{n}}^{2}\\ \vdots & \vdots & \ddots & \vdots \\ S_{t_{0}}^{d} & S_{t_{1}}^{d} & \cdots & S_{t_{n}}^{d}\\ \end{pmatrix} = \begin{pmatrix} S_0 & S_{t_{1}}^{1} & \cdots & S_{t_{n}}^{1}\\ S_0 & S_{t_{1}}^{2} & \cdots & S_{t_{n}}^{2}\\ \vdots & \vdots & \ddots & \vdots \\ S_0 & S_{t_{1}}^{d} & \cdots & S_{t_{n}}^{d}\\ \end{pmatrix} $$
St = prixActif(t,W,sigma_input,r_input,S0_input)
plt.subplots(figsize=(8,5))
for i in range(4):
plt.plot(t,St[i], label = "actif n°{}".format(i+1))
plt.xlabel("temps")
plt.ylabel("prix")
plt.title("Evolution du prix d'actifs [BS] | $\sigma$ = {}, n = {}".format(sigma_input,n_input))
plt.legend()
plt.show()
def asian_payoff(Mat_st,K = K_input):
"""
Fonction permettant d'évaluer le payoff de l'option asiatique.
--------
Arguments:
Mat_st : Matrice des trajectoires de l'actif.
K : Strike.
Returns:
payoff : Vecteur des payoffs pour chacune des trajectoires de l'actif (i.e chacune des lignes de Mat_st par construction que W).
"""
payoff = np.mean(Mat_st, axis = 1) - K
payoff = payoff * (payoff >= 0)
return payoff
payoff_temp = asian_payoff(St,K_input)
for i in range(4):
print("payoff = {}".format(payoff_temp[i]))
payoff = 7.387975617734099 payoff = 1.1241555517631099 payoff = 46.0255750009068 payoff = 22.897497405783625
Jusque-là on a seulement calculer une réalisation du payoff de notre indice. Toutefois on veut le prix de l'option, ce qui correspond à $ \mathbb{E} \Bigl[ \bigl(I_{T} - K \bigr)_{+} \Bigr] $.
Pour cela on va passer par la méthode de Monte-Carlo. Cette méthode est basée sur la loi forte des grands nombres car on sait que pour toute variable intégrable $Z$, si on a des variables $Z_{i}$ de même loi que $Z$, on a $\overline{Z}_{M}= \frac{1}{M} \sum \limits_{\underset{}{i=0}}^M Z_i$ $ \underset{M \to +\infty}{\overset{p.s}{\longrightarrow}} \mathbb{E}[Z]$.
Cette méthode étant probabiliste, il nous fait un intervalle de confiance de notre estimateur de $\mathbb{E}[Z]$.
Par le TCL (Théorème Central Limite), on sait que (en rajoutant l'hypothèse que Z soit carré intégrable) $\sqrt{M} \frac{ \overline{Z}_{M}-\mathbb{E}[Z] }{ \sqrt{\mathbb{V}(Z)} } \sim \mathcal{N}(0,1)$. De ce fait on peut construire l'intervalle suivant pour $\mathbb{E}[Z]$ :
$$
\begin{equation*}
\Bigl[ \overline{Z}_{M} - q_{\alpha}\sqrt{\frac{\mathbb{V}(Z)}{M}} ; \overline{Z}_{M} + q_{\alpha} \sqrt{\frac{\mathbb{V}(Z)}{M}} \Bigr] \\
\text{avec } q_{\alpha} \text{le quantile d'ordre } \alpha \text{ de la loi normale centrée réduite}
\end{equation*}
$$
Ici $Z = \bigl(I_{T} - K \bigr)_{+}$ et comme on ne connaît pas la loi du payoff, on déterminera $\mathbb{V} \Bigl[ \bigl(I_{T} - K \bigr)_{+} \Bigr]$ numériquement à partir des tirages déjà fais de $\bigl(I_{T} - K \bigr)_{+}$.
def IntervalleConfiance(tirages, quantile = 1.96):
"""
Fonction calculant un intervalle de confiance à 95% pour la moyenne d'un vecteur tirages.
--------
Arguments:
tirages : Vecteur des données dont la moyenne est à estimer.
quantile : Quantile de la loi Normale centrée réduite pour un niveau de confiance donné. (Par défaut 1.96 (95%))
Returns:
I : Intervalle d'estimation de la moyenne pour un niveau de confiance donné. (Par défaut 1.96 (95%))
"""
I = np.zeros(2)
I += tirages.mean()
var = tirages.var()
Mc = len(tirages)
I[0] = I[0] - quantile*np.sqrt(var/Mc)
I[1] = I[1] + quantile*np.sqrt(var/Mc)
return(I)
def prix_option_asiatique(S0 = S0_input, K = K_input, r = r_input, sigma = sigma_input, T = T_input, n = n_input, Mc = Mc_input):
"""
Fonction estimant le prix d'une option asiatique à partir de la méthode de Monte-Carlo.
--------
Arguments:
S0 : Spot.
K : Strike.
r : Taux d'intérêt de l'actif sans risque.
sigma : Volatilité de l'actif.
T : Maturité.
n : Nombre de points de la discrétisation.
Mc : Nombre d'échantillons pour la méthode de Monte-Carlo.
Returns:
prix : Prix de l'option estimé.
intervalle : Intervalle de confiance pour l'estimation du prix. (seuil à 95% par défaut).
"""
t,Mb = dMouvementBrownien(T,n,Mc)
St = prixActif(t,Mb,sigma,r,S0)
payoff_vec = asian_payoff(St,K)
intervalle = np.exp(-r*T)*IntervalleConfiance(payoff_vec)
prix = np.exp(-r*T)*payoff_vec.mean()
return(prix,intervalle)
t1 = time.time()
prix,intervalle = prix_option_asiatique()
t2 = time.time()
print("Avec n = {} et Mc = {}".format(n_input,Mc_input))
print("On trouve les valeurs suivantes :")
print("borne inf intervalle 95% : {}".format(round(intervalle[0],3)))
print("prix : {}".format(round(prix,3)))
print("borne sup intervalle 95% : {}".format(round(intervalle[1],3)))
print("temps de calcul : {}s".format(round(t2-t1,1)))
Avec n = 200 et Mc = 10000 On trouve les valeurs suivantes : borne inf intervalle 95% : 8.633 prix : 8.816 borne sup intervalle 95% : 8.999 temps de calcul : 0.2s
MonteCarlo = np.array([100,200,500,800,1000,5000,10000,50000,100000,500000,1000000])
prix = np.zeros(len(MonteCarlo))
Intervalle = np.zeros( (len(MonteCarlo),2) )
temps = []
print("Pour {} points de discrétisation".format(n_input+1))
for i in range(len(MonteCarlo)):
t1 = time.time()
prix[i],Intervalle[i] = prix_option_asiatique(Mc = MonteCarlo[i])
t2 = time.time()
temps.append(t2-t1)
if( (t2-t1)/60 < 1 ):
print("pour {} itérations de Monte-Carlo, temps de calcul : {}s".format( MonteCarlo[i], round(t2-t1,2) ))
else :
print("pour {} itérations de Monte-Carlo, temps de calcul : {}min".format( MonteCarlo[i], round( (t2-t1)/60,2) ) )
Pour 201 points de discrétisation pour 100 itérations de Monte-Carlo, temps de calcul : 0.01s pour 200 itérations de Monte-Carlo, temps de calcul : 0.04s pour 500 itérations de Monte-Carlo, temps de calcul : 0.01s pour 800 itérations de Monte-Carlo, temps de calcul : 0.02s pour 1000 itérations de Monte-Carlo, temps de calcul : 0.02s pour 5000 itérations de Monte-Carlo, temps de calcul : 0.09s pour 10000 itérations de Monte-Carlo, temps de calcul : 0.18s pour 50000 itérations de Monte-Carlo, temps de calcul : 0.87s pour 100000 itérations de Monte-Carlo, temps de calcul : 1.83s pour 500000 itérations de Monte-Carlo, temps de calcul : 9.84s pour 1000000 itérations de Monte-Carlo, temps de calcul : 22.9s
fig = plt.figure(figsize=(16,5))
ax1 = fig.add_subplot(121); ax2 = fig.add_subplot(122)
ax1.scatter(MonteCarlo,prix, color = "k", label = "Prix")
# ax1.plot(MonteCarlo, Intervalle[:,0], color = 'r', label = "intervalle de confiance 95%")
# ax1.plot(MonteCarlo, Intervalle[:,1], color = 'r')
ax1.fill_between(MonteCarlo, Intervalle[:,0], Intervalle[:,1], alpha=.2, linewidth=0, color='grey', label = "intervalle de confiance 95%")
ax1.set_title("Evolution du prix de l'option asiatique | $\sigma$ = {}, K = {}".format(sigma_input,K_input))
ax1.set_xlabel("Nombre d'itérations de MonteCarlo (échelle log)")
ax1.set_ylabel("Prix de l'option")
ax1.set_xscale('log')
ax1.set_xlim(100,1000000)
ax1.set_ylim(4,12)
ax1.grid()
ax2.plot(MonteCarlo, temps, label = "Temps de calcul", color='k')
ax2.set_ylabel("Temps de calcul (secondes)")
ax2.set_xlabel("Nombre d'itérations de Monte-Carlo")
ax2.set_title("Evolution du temps de calcul")
ax2.set_xscale('log')
ax2.set_xlim(100,1000000)
ax2.set_ylim(0,25)
ax2.grid()
ax1.legend(); ax2.legend(); plt.show()
$\mathbf{Q1)}$
$\mathbf{Q2})$
## variable de contrôle non actualisé (i.e e^{r*T} Y )
def control_variate1(Mat_St):
return(np.exp(np.mean(np.log(Mat_St), axis = 1)))
def asian_option_controlled(S0 = S0_input, K = K_input, r = r_input, sigma = sigma_input, T = T_input, n = n_input, Mc = Mc_input):
"""
Fonction estimant le prix d'une option asiatique à partir de la méthode de Monte-Carlo et à l'aide d'une variable de contrôle.
Ici la variable de contrôle est :
Y = exp^{-rT}exp{\frac{1}{T}\int_0^T log(S_t) dt}
--------
Arguments:
S0 : Spot.
K : Strike.
r : Taux d'intérêt de l'actif sans risque.
sigma : Volatilité de l'actif.
T : Maturité.
n : Nombre de points de la discrétisation.
Mc : Nombre d'échantillons pour la méthode de Monte-Carlo.
Returns:
prix : Prix de l'option estimé.
intervalle : Intervalle de confiance pour l'estimation du prix. (seuil à 95% par défaut).
"""
t,Mb = dMouvementBrownien(T,n,Mc)
St = prixActif(t,Mb,sigma,r,S0)
payoff_vec = asian_payoff(St,K)
Y = control_variate1(St)
Mat_cov = np.cov(payoff_vec,Y)
Beta_opti = Mat_cov[0,1]/Mat_cov[1,1]
esp_Y = S0*np.exp(- T/2 * (r + (sigma**2)/ 6) )
Q = payoff_vec - Beta_opti*Y
intervalle = np.exp(-r*T)*IntervalleConfiance(Q) + Beta_opti*esp_Y
prix = np.exp(-r*T)*Q.mean() + Beta_opti*esp_Y
return(prix,intervalle)
t1 = time.time()
p,I = asian_option_controlled()
t2 = time.time()
print("Avec n = {} et Mc = {}".format(n_input,Mc_input))
print("On trouve les valeurs suivantes :")
print("borne inf intervalle 95% : {}".format(round(I[0],3)))
print("prix : {}".format(round(p,3)))
print("borne sup intervalle 95% : {}".format(round(I[1],3)))
print("temps de calcul : {}s".format(round(t2-t1,1)))
Avec n = 200 et Mc = 10000 On trouve les valeurs suivantes : borne inf intervalle 95% : 8.718 prix : 8.772 borne sup intervalle 95% : 8.826 temps de calcul : 0.2s
def Comparaison_naif_control_asiatique(S0 = S0_input, K = K_input, r = r_input, sigma = sigma_input, T = T_input, n = n_input, Mc = Mc_input):
""" Fonction comparant les estimations des prix pour deux méthodes différentes à partir d'un même processus. """
t,Mb = dMouvementBrownien(T,n,Mc)
St = prixActif(t,Mb,sigma,r,S0)
payoff_vec = asian_payoff(St,K)
Y = control_variate1(St)
Mat_cov = np.cov(payoff_vec,Y)
Beta_opti = Mat_cov[0,1]/Mat_cov[1,1]
esp_Y = S0*np.exp(- T/2 * (r + (sigma**2)/ 6) )
Q = payoff_vec - Beta_opti*Y
intervalle_naif = np.exp(-r*T)*IntervalleConfiance(payoff_vec)
prix_naif = np.exp(-r*T)*payoff_vec.mean()
intervalle_control = np.exp(-r*T)*IntervalleConfiance(Q) + Beta_opti*esp_Y
prix_control = np.exp(-r*T)*Q.mean() + Beta_opti*esp_Y
return(prix_naif,intervalle_naif,prix_control,intervalle_control)
vol_sigma = np.linspace(0.2,2,10)
prix_naif = np.zeros(len(vol_sigma)) ; prix_control = np.zeros(len(vol_sigma))
Intervalle_naif = np.zeros( (len(vol_sigma),2) ) ; Intervalle_control = np.zeros( (len(vol_sigma),2) );
temps = 0
print("Pour {} points de discrétisation et Mc = {}".format(n_input+1,Mc_input))
t1 = time.time()
for i in range(len(vol_sigma)):
prix_naif[i],Intervalle_naif[i],prix_control[i],Intervalle_control[i] = Comparaison_naif_control_asiatique(sigma = vol_sigma[i])
t2 = time.time()
print("temps de calcul : {}s".format(round(t2-t1,3)))
longueur_naif = np.zeros(len(vol_sigma))
longueur_control = np.zeros(len(vol_sigma))
for i in range(len(vol_sigma)):
longueur_naif[i] = Intervalle_naif[i,1]-Intervalle_naif[i,0]
longueur_control[i] = Intervalle_control[i,1]-Intervalle_control[i,0]
Pour 201 points de discrétisation et Mc = 10000 temps de calcul : 1.688s
fig = plt.figure(figsize=(16,5))
ax1 = fig.add_subplot(121); ax2 = fig.add_subplot(122)
# ax1.scatter(vol_sigma,prix_naif, color = "b", label = "prix naif")
# ax1.plot(vol_sigma, Intervalle_naif[:,0], color = 'r', label = "intervalle de confiance 95% [naif]")
# ax1.plot(vol_sigma, Intervalle_naif[:,1], color = 'r')
# ax1.scatter(vol_sigma,prix_control, color = "c", label = "prix control")
# ax1.plot(vol_sigma, Intervalle_control[:,0], color = 'm', label = "intervalle de confiance 95% [control]")
# ax1.plot(vol_sigma, Intervalle_control[:,1], color = 'm')
ax1.errorbar(vol_sigma, prix_naif, longueur_naif/2, fmt='o', color='k', ecolor='black', linestyle='dotted', linewidth=1, capsize=3, label="Prix naif (bornes à 95%)")
ax1.errorbar(vol_sigma, prix_control, longueur_control/2, fmt='v', color='mediumslateblue', ecolor='mediumslateblue', linestyle='dashed', linewidth=1, capsize=3, label="Prix contrôlé par Y (bornes à 95%)")
ax1.fill_between(vol_sigma, np.minimum(Intervalle_naif[:,0],Intervalle_control[:,0]), np.maximum(Intervalle_naif[:,1],Intervalle_control[:,1]), alpha=.2, linewidth=0, color='grey')
ax1.set_title("Evolution du prix de l'option asiatique | Mc = {}".format(Mc_input))
ax1.set_xlabel("Volatilité de l'actif risqué")
ax1.set_ylabel("Prix de l'option")
ax1.grid()
ax2.plot(vol_sigma,longueur_naif, label = "Taille intervalle naïf", color='k', linewidth=1)
ax2.plot(vol_sigma,longueur_control, label = "Taille intervalle [contrôle avec $Y$]", color='mediumslateblue', linestyle='dashed', linewidth=1)
ax2.set_title("Comparaison des tailles des intervalles de confiance à $95\%$ | Mc = {}".format(Mc_input));
ax2.set_xlabel("Volatilité de l'actif risqué")
ax2.set_ylabel("Longueur")
ax2.set_ylim(0,16)
ax2.grid()
ax1.legend(); ax2.legend(); plt.show()
plt.show()
$\mathbf{Q3})$
## variable de contrôle non actualisé (i.e e^{r*T} Z )
def control_variate2(Mat_St,K = K_input):
control_variate = np.exp(np.mean(np.log(Mat_St), axis = 1)) - K
control_variate = control_variate * (control_variate >= 0)
return control_variate
def asian_option_controlled_2(S0 = S0_input, K = K_input, r = r_input, sigma = sigma_input, T = T_input, n = n_input, Mc = Mc_input):
"""
Fonction estimant le prix d'une option asiatique à partir de la méthode de Monte-Carlo et à l'aide d'une variable de contrôle.
Ici la variable de contrôle est :
Z = exp^{-rT}(exp{\frac{1}{T}\int_0^T log(S_t) dt} - K)_+
--------
Arguments:
S0 : Spot.
K : Strike.
r : Taux d'intérêt de l'actif sans risque.
sigma : Volatilité de l'actif.
T : Maturité.
n : Nombre de points de la discrétisation.
Mc : Nombre d'échantillons pour la méthode de Monte-Carlo.
Returns:
prix : Prix de l'option estimé.
intervalle : Intervalle de confiance pour l'estimation du prix. (seuil à 95% par défaut).
"""
t,Mb = dMouvementBrownien(T,n,Mc)
St = prixActif(t,Mb,sigma,r,S0)
payoff_vec = asian_payoff(St,K)
Z = control_variate2(St,K)
Mat_cov = np.cov(payoff_vec,Z)
Beta_opti = Mat_cov[0,1]/Mat_cov[1,1]
X = np.random.normal(size = int(2*Mc))
d = 1/sigma * np.sqrt(3/T) * (np.log(K/S0) - (r- sigma**2/2)*T/2)
N1 = (X <= sigma*np.sqrt(T/3) - d).mean()
N2 = (X <= -d).mean()
esp_Z = S0*np.exp(- T/2 * (r + (sigma**2)/6) )*N1 - K*np.exp(-r*T)*N2
Q = payoff_vec - Beta_opti*Z
intervalle = np.exp(-r*T)*IntervalleConfiance(Q) + Beta_opti*esp_Z
prix = np.exp(-r*T)*Q.mean() + Beta_opti*esp_Z
return(prix,intervalle)
t1 = time.time()
p,I = asian_option_controlled_2()
t2 = time.time()
print("Avec n = {} et Mc = {}".format(n_input,Mc_input))
print("On trouve les valeurs suivantes :")
print("borne inf intervalle 95% : {}".format(round(I[0],3)))
print("prix : {}".format(round(p,3)))
print("borne sup intervalle 95% : {}".format(round(I[1],3)))
print("temps de calcul : {}s".format(round(t2-t1,1)))
Avec n = 200 et Mc = 10000 On trouve les valeurs suivantes : borne inf intervalle 95% : 8.874 prix : 8.878 borne sup intervalle 95% : 8.883 temps de calcul : 0.2s
Comparatif des variables de contrôle
def Comparaison_asiatique(S0 = S0_input, K = K_input, r = r_input, sigma = sigma_input, T = T_input, n = n_input, Mc = Mc_input):
""" Fonction de comparaison des estimations des prix à partir d'un même processus. """
t,Mb = dMouvementBrownien(T,n,Mc)
St = prixActif(t,Mb,sigma,r,S0)
payoff_vec = asian_payoff(St,K)
Y = control_variate1(St)
Z = control_variate2(St,K)
Mat_cov1 = np.cov(payoff_vec,Y)
Mat_cov2 = np.cov(payoff_vec,Z)
Beta_opti_Y = Mat_cov1[0,1]/Mat_cov1[1,1]
Beta_opti_Z = Mat_cov2[0,1]/Mat_cov2[1,1]
X = np.random.normal(size = int(2*Mc))
d = 1/sigma * np.sqrt(3/T) * (np.log(K/S0) - (r- sigma**2/2)*T/2)
N1 = (X <= sigma*np.sqrt(T/3) - d).mean()
N2 = (X <= -d).mean()
esp_Z = S0*np.exp(- T/2 * (r + (sigma**2)/6) )*N1 - K*np.exp(-r*T)*N2
esp_Y = S0*np.exp(- T/2 * (r + (sigma**2)/ 6) )
Q1 = payoff_vec - Beta_opti_Y*Y
Q2 = payoff_vec - Beta_opti_Z*Z
intervalle_naif = np.exp(-r*T)*IntervalleConfiance(payoff_vec)
prix_naif = np.exp(-r*T)*payoff_vec.mean()
intervalle_control1 = np.exp(-r*T)*IntervalleConfiance(Q1) + Beta_opti_Y*esp_Y
prix_control1 = np.exp(-r*T)*Q1.mean() + Beta_opti_Y*esp_Y
intervalle_control2 = np.exp(-r*T)*IntervalleConfiance(Q2) + Beta_opti_Z*esp_Z
prix_control2 = np.exp(-r*T)*Q2.mean() + Beta_opti_Z*esp_Z
return(prix_naif,intervalle_naif,prix_control1,intervalle_control1,prix_control2,intervalle_control2)
vol_sigma = np.linspace(0.2,2,10)
prix_naif = np.zeros(len(vol_sigma)) ; prix_control1 = np.zeros(len(vol_sigma)) ; prix_control2 = np.zeros(len(vol_sigma));
Intervalle_naif = np.zeros( (len(vol_sigma),2) ) ; Intervalle_control1 = np.zeros( (len(vol_sigma),2) ); Intervalle_control2 = np.zeros( (len(vol_sigma),2) );
temps = 0
print("Pour {} points de discrétisation et Mc = {}".format(n_input+1,Mc_input))
t1 = time.time()
for i in range(len(vol_sigma)):
prix_naif[i],Intervalle_naif[i],prix_control1[i],Intervalle_control1[i],prix_control2[i],Intervalle_control2[i] = Comparaison_asiatique(sigma = vol_sigma[i])
t2 = time.time()
print("temps de calcul : {}s".format(round(t2-t1,3)))
longueur_naif = np.zeros(len(vol_sigma))
longueur_control1 = np.zeros(len(vol_sigma))
longueur_control2 = np.zeros(len(vol_sigma))
for i in range(len(vol_sigma)):
longueur_naif[i] = Intervalle_naif[i,1]-Intervalle_naif[i,0]
longueur_control1[i] = Intervalle_control1[i,1]-Intervalle_control1[i,0]
longueur_control2[i] = Intervalle_control2[i,1]-Intervalle_control2[i,0]
Pour 201 points de discrétisation et Mc = 10000 temps de calcul : 1.739s
fig = plt.figure(figsize=(16,5))
ax1 = fig.add_subplot(121); ax2 = fig.add_subplot(122)
# ax1.scatter(vol_sigma,prix_naif, color = "b", label = "prix naif")
# ax1.plot(vol_sigma, Intervalle_naif[:,0], color = 'r', label = "intervalle de confiance 95% [naif]")
# ax1.plot(vol_sigma, Intervalle_naif[:,1], color = 'r')
ax1.errorbar(vol_sigma, prix_naif, longueur_naif/2, fmt='o', color='k', ecolor='black', linestyle='dotted', linewidth=1, capsize=3, label="Prix naif (bornes à 95%)")
# ax1.scatter(vol_sigma,prix_control1, color = "c", label = "prix control - Y")
# ax1.plot(vol_sigma, Intervalle_control1[:,0], color = 'm', label = "intervalle de confiance 95% [control - Y]")
# ax1.plot(vol_sigma, Intervalle_control1[:,1], color = 'm')
ax1.errorbar(vol_sigma, prix_control1, longueur_control1/2, fmt='v', color='mediumslateblue', ecolor='mediumslateblue', linestyle='dashed', linewidth=1, capsize=3, label="Prix contrôlé par Y (bornes à 95%)")
# ax1.scatter(vol_sigma,prix_control2, color = "y", label = "prix control - Z")
# ax1.plot(vol_sigma, Intervalle_control2[:,0], color = 'g', label = "intervalle de confiance 95% [control - Z]")
# ax1.plot(vol_sigma, Intervalle_control2[:,1], color = 'g')
ax1.errorbar(vol_sigma, prix_control2, longueur_control2/2, fmt='d', color='royalblue', ecolor='royalblue', linestyle='solid', linewidth=1, capsize=3, label="Prix contrôlé par Z (bornes à 95%)")
ax1.set_title("Evolution du prix de l'option asiatique | Mc = {}".format(Mc_input))
ax1.set_xlabel("Volatilité de l'actif risqué")
ax1.set_ylabel("Prix de l'option")
ax1.grid()
ax2.plot(vol_sigma,longueur_naif, label = "Taille intervalle [naïf]", color='k', linestyle='dotted')
ax2.plot(vol_sigma,longueur_control1, label = "Taille intervalle [contrôle avec - $Y$]", color='mediumslateblue', linestyle='dashed')
ax2.plot(vol_sigma,longueur_control2, label = "Taille intervalle [contrôle avec - $Z$]", color='royalblue', linestyle='solid')
ax2.set_title("Comparaison des tailles des intervalles de confiance à $95\%$ | Mc = {}".format(Mc_input));
ax2.set_xlabel("Volatilité de l'actif risqué")
ax2.set_ylabel("Longueur")
ax2.set_ylim(0,16)
ax2.grid()
ax1.legend(); ax2.legend(); plt.show()
plt.show()
Nos techniques de réduction de varance marchent bien. Les graphes ci-dessus nous montrent bien que les méthodes permettent de converger vers le même prix (donc c'est bien implémenté) et aussi les méthodes de réduction de variance réduisent effectivement la variance de nos estimateurs car la longueur de l'intervalle de confiance se réduit avec les variables de contrôle et de plus l'utilisation de $Z$ nous donne de meilleurs résultats qu'avec $Y$.
$\mathbf{Q5})$
def Compute_VegaPWMethod(S, W, r, sg, Tdiscret, K):
"""
Fonction calculant le Vega d'une option asiatique par la méthode pathwise. Le prix P de l'option est défini comme
P = e^{-rT} E[ (\bar{S} - K)_+ ] , où \bar{S} désigne la moyenne de S.
--------
Arguments :
S : Sous-jacent S, matrice contenant M trajectoires de S. M désignant le nombre de tirages de Monte-Carlo.
W : Mouvement Brownien associé à S. En particulier np.shape(S)=np.shape(W).
r : Taux d'intérêt de l'actif sans-risque.
sg : Volatilité du sous-jacent.
Tdiscret : Discrétisation spatiale de l'option.
K : Strike.
Return :
V : Vega de l'option asiatique.
"""
MC = np.shape(S)[0]
#On cherche les trajectoires qui nous intéressent
S_sel = np.mean(S, axis=1)>K #moyenne selon l'axe 1 si les trajectoires sont en ligne.
#On calcule l'espérance
Tscal = np.outer(np.ones(MC),Tdiscret) #On créer une matrice de tps: on crée autant de ligne de dsicrétisation temporelle qu'il y a de trajectoires.
Mderiv = np.mean((W-sg*Tscal)*S, axis=1) #On fait la moyenne qui se trouve dans l'espérance
#Vérification
assert(np.shape(Mderiv)==np.shape(S_sel))
assert(np.size(S_sel)==MC)
#On annule les valeurs ne qui ne nous intéresse pas
for i in range(MC):
if (S_sel[i]==0):
Mderiv[i]=0
#Monte-Carlo
#Vega
Vtirages = np.exp(-r*Tdiscret[-1])*Mderiv
V = np.mean(Vtirages)
return V, Vtirages
def Compute_VegaLLHMethod(S, r, sg, Tdiscret, K):
"""
Fonction calculant le Vega d'une option asiatique par la méthode pathwise. Le prix P de l'option est défini comme
P = e^{-rT} E[ (\bar{S} - K)_+ ] , où \bar{S} désigne la moyenne de S.
--------
Arguments :
S : Sous-jacent S, matrice contenant M trajectoires de S. M désignant le nombre de tirages de Monte-Carlo.
W : Mouvement Brownien associé à S. En particulier np.shape(S)=np.shape(W).
r : Taux d'intérêt de l'actif sans-risque.
sg : Volatilité du sous-jacent.
Tdiscret : Discrétisation spatiale de l'option.
K : Strike.
Return :
V : Vega de l'option asiatique.
"""
MC = np.shape(S)[0]
#On cherche les trajectoires qui nous intéressent
S_sel = np.mean(S, axis=1)>K #moyenne selon l'axe 1 si les trajectoires sont en ligne.
#On calcul le terme somme dans l'espérance, en simultané pour chaque trajectoire
Svega = np.zeros(MC)
for i in range(np.shape(S)[1]-1):
# Préliminaires
zeta = (np.log(S[:,i+1]/S[:,i]) - (r - (sg**2)/2) * (Tdiscret[i+1] - Tdiscret[i])) / (sg * np.sqrt(Tdiscret[i+1] - Tdiscret[i]))
# Vérification
assert(np.size(zeta)==np.size(Svega))
# Calcul somme
Svega += -zeta *(np.sqrt(Tdiscret[i+1] - Tdiscret[i]) - zeta/sg) - 1/sg
#On annule les valeurs ne qui ne nous intéresse pas
for i in range(MC):
if (S_sel[i]==0):
Svega[i]=0
#Monte-Carlo
# E = np.mean(Svega)
#Vega
Vtirages = np.exp(-r*Tdiscret[-1])*Svega
V = np.mean(Vtirages)
return V, Vtirages
##### A valider
t,Mb = dMouvementBrownien(T_input,n_input,Mc_input)
St = prixActif(t,Mb,sigma_input,r_input,S0_input)
VegaPW, VegaPWtirages = Compute_VegaPWMethod(St, Mb, r_input, sigma_input, t, K_input)
VegaLLH, VegaLLHtirages = Compute_VegaLLHMethod(St, r_input, sigma_input, t, K_input)
IC_VegaPW = IntervalleConfiance(VegaPWtirages)
IC_VegaLLH = IntervalleConfiance(VegaLLHtirages)
print('Vega Pathwise : ', VegaPW)
print('Vega Pathwise IC95% : ', IC_VegaPW)
print('Vega LLH : ', VegaLLH)
print('Vega LLH IC95% : ', IC_VegaLLH)
Vega Pathwise : 17.707533732137144 Vega Pathwise IC95% : [16.84911659 18.56595087] Vega LLH : -1.5453554297335932 Vega LLH IC95% : [-3.14395058 0.05323972]
#Paramètres
S0 = 100
sg = 0.3
r = 0.05
T = 4
t1 = T/4
coef = 0.9
MC = 30000
$\mathbf{Q1})$
avec $d_1 = \frac{log(\frac{1}{\lambda})+(r+\frac{\sigma^2}{2})(T-t_1)}{\sigma \sqrt{T-t_1}}$ et $d_2 = d_1 - \sigma \sqrt{T-t_1}$.
$\mathbf{Q2})$
def MC4Computing(S0, sg, r, T, t1, coef, MC):
"""
Fonction estimant le prix d'une delayed start option, à partir de l'algorithme 1,
donné dans le sujet.
--------
Arguments :
S0 : Spot de l'actif.
sg : Volatilité de l'actif.
r : Rendement de l'actif sans risque.
T : Maturité.
t1 : Temps d'intérêt pour la delayed start option.
coef : Coefficient d'intérêt pour la delayed start option.
MC : Nombre d'itérations de Monte-Carlo
Returns :
DP : Estimation par Monte-Carlo du prix de la delayed start option.
payoffs : Les payoffs actualisés ayant servi à l'estimation.
"""
#Initialisation
payoffs = np.zeros(MC)
for i in range(MC):
#1
Z1 = np.random.normal()
St1 = S0 * np.exp((r - (sg**2) /2) * t1 + sg * np.sqrt(t1) * Z1)
#2
K = coef * St1
#3
Z2 = np.random.normal()
ST = St1 * np.exp((r - (sg**2) /2) * (T - t1) + sg * np.sqrt(T - t1) * Z2)
#4
payoffs[i] = np.exp(-r * T) * max(ST - K, 0)
DP = np.mean(payoffs)
return DP, payoffs
def analytical_F(S,K,sg,r,T,t):
"""
Calcul de F(St, K, sg, r, T-t) = S_t N(d_1) - K exp(-r (T-t)) N(d_2)
"""
d1 = (np.log(S/K) + (r + sg**2/2) * (T-t)) / (sg * np.sqrt(T-t))
d2 = d1 - sg * np.sqrt(T-t)
return S * norm.cdf(d1) - K * np.exp(-r * (T-t)) * norm.cdf(d2)
def Conditional_MC4Computing(S0, sg, r, T, t1, coef, MC):
"""
Fonction estimant le prix d'une delayed start option, à partir de l'algorithme 2,
donné dans le sujet.
--------
Arguments :
S0 : Spot de l'actif.
sg : Volatilité de l'actif.
r : Rendement de l'actif sans risque.
T : Maturité.
t1 : Temps d'intérêt pour la delayed start option.
coef : Coefficient d'intérêt pour la delayed start option.
MC : Nombre d'itérations de Monte-Carlo
Returns :
DP : Estimation par Monte-Carlo (conditionnel) du prix de la delayed start option.
payoffs : Les payoffs actualisés ayant servi à l'estimation.
"""
#Initialisation
payoffs = np.zeros(MC)
for i in range(MC):
#1
Z1 = np.random.normal()
St1 = S0 * np.exp((r - (sg**2) /2) * t1 + sg * np.sqrt(t1) * Z1)
#2
K = coef * St1
#3
payoffs[i] = np.exp(-r * t1) * analytical_F(St1, K, sg, r, T, t1)
DP = np.mean(payoffs)
return DP, payoffs
DP_classique, Pclassique = MC4Computing(S0, sg, r, T, t1, coef, MC)
DP_cond, Pcond = Conditional_MC4Computing(S0, sg, r, T, t1, coef, MC)
print("Calssique : ", DP_classique)
print("Cond : ", DP_cond)
Calssique : 31.889599614105208 Cond : 31.447204167665216
vol_sigma = np.linspace(0.2,2.,10)
DP_classique = np.zeros(10)
LIclassique = np.zeros(10)
DP_cond = np.zeros(10)
LIcond = np.zeros(10)
for i in range(len(vol_sigma)):
#Compute PRICE
DP_classique[i], Pclassique = MC4Computing(S0, vol_sigma[i], r, T, t1, coef, MC)
DP_cond[i], Pcond = Conditional_MC4Computing(S0, vol_sigma[i], r, T, t1, coef, MC)
#Compute Intervalle
Iclassique = IntervalleConfiance(Pclassique)
Icond = IntervalleConfiance(Pcond)
#Compute Longueur
LIclassique[i] = Iclassique[1] - Iclassique[0]
LIcond[i] = Icond[1] - Icond[0]
fig = plt.figure(figsize=(16,5))
ax1 = fig.add_subplot(121); ax2 = fig.add_subplot(122)
ax1.errorbar(vol_sigma, DP_classique, LIclassique/2, fmt='o', color='royalblue', ecolor='royalblue', linestyle='dotted', linewidth=1, capsize=3, label="Prix classique (bornes à 95%)")
ax1.errorbar(vol_sigma, DP_cond, LIcond/2, fmt='v', color='mediumslateblue', ecolor='mediumslateblue', linestyle='dashed', linewidth=1, capsize=3, label="Prix conditionné (bornes à 95%)")
ax1.set_title("Evolution du prix de la delayed start option | Mc = {}".format(MC))
ax1.set_xlabel("Volatilité de l'actif risqué")
ax1.set_ylabel("Prix de l'option")
ax1.grid()
ax2.plot(vol_sigma, LIclassique, label = "Taille intervalle [classique]", color='royalblue', linestyle='dotted')
ax2.plot(vol_sigma, LIcond, label = "Taille intervalle [conditionné]", color='mediumslateblue', linestyle='dashed')
ax2.set_title("Comparaison des tailles des intervalles de confiance à $95\%$ | Mc = {}".format(MC));
ax2.set_xlabel("Volatilité de l'actif risqué")
ax2.set_ylabel("Longueur")
ax2.grid()
ax1.legend(); ax2.legend(); plt.show()
plt.show()
Ici notre objectif est de déterminer le prix d'un put bermudien : $$ \begin{equation*} p = \underset{\tau \in \{t_0,..,t_n \} }{sup} \mathbb{E} \left[ e^{-r \tau} \left( K - S_{\tau} \right)^{+} \right] \\ \end{equation*} $$
# Paramètres
S0 = 90
K = 100
r = 0.03
sg = 0.15
T = 1
# Discretisation
n = 200
# Monte-Carlo
MC = 1000
# Payoff de l'option étudiée
f_payoff = lambda S : np.maximum(K-S,0)
$\mathbf{Q1})$
def Compute_H(S,Basis):
"""
Fonction calculant la matrice H de l'algorithme de Longstaff-Schwarz et définie comme
\forall l,k \in [1,N]^2, H_{l,k} = \sum_{j=1}^{MC} Basis_l(S_j) * Basis_k(S_j)
où S est le vecteur des j trajectoires de Monte-Carlo de l'actif S prisent au temps k fixé.
La matrice H est donc dépendante du temps k.
--------
Arguments:
S : Vecteurs des valeurs à un instant donné des trajectoires de S. Taille : (MC,1) où MC est le nombre d'itérations de Monte-Carlo.
Basis : Base fonctionnelle de notre algorithme de Longstaff-Schwarz (Numpy Array Python).
Returns:
H : Matrice au temps k utile pour l'algorithme de Longstaff-Schwarz.
"""
# Déclaration et Initialisation
L = len(Basis)
MC = np.size(S)
H = np.zeros((L,L))
P = np.zeros((MC,L))
# Remplissage de P
for i in range(L):
P[:,i] = Basis[i](S)
P = np.array(P)
# Remplissage de H
for k in range(L):
Pk = np.outer(P[:,k],np.ones(L))
Hrow_k = P * Pk
H[k,:] = np.sum(Hrow_k, axis=0) #H_{k,l} = H_{l,k} donc H[k,:] = H[:,k]
return H
def Compute_alpha(S,S_tau_star,Basis,f_payoff,t,tau_star,r):
"""
Fonction calculant le vecteur \alpha de L coordonnées à l'instant t d'après la formule matricielle
H \alpha = Pbis
avec H la matrice à l'instant t obtenue par la fonction Compute_H, et Pbis le vecteur de L coordonnées défini comme
\forall l \in [1,L], (Pbis)_{l} = \sum_{j=1}^{MC} e^{-r * (\tau_star_{j} - t)} * f_payoff(S_{\tau_star_{j}}) * P_{l}(S_{t}^{j}).
On note P_{l}, le l-ième élément de notre base fonctionnelle Basis.
--------
Arguments:
S : Vecteurs des valeurs, à un instant t donné, des trajectoires de S. Taille : (MC,1) où MC est le nombre d'itérations de Monte-Carlo.
S_tau_star : Valeurs, à un instant tau_star, de la j-ieme trajectoires de S.
Basis : Base fonctionnelle de notre algorithme de Longstaff-Schwarz (Numpy Array Python).
f_payoff : Payoff de l'option étudiée.
t : Instant auquel on souhaite obtenir le vecteur \alpha minimisant une certaine distance.
tau_star : Valeurs de tau_star à l'instant t+1.
r : Rendement de l'actif sans-risque.
Returns:
\alpha : Vecteur des \alpha à l'instant t. Taille : (L,1) où L est le nombre d'éléments de la base fonctionnelle Basis.
"""
#Déclaration et Initialisation
L = len(Basis)
n = np.size(S)
Pbis = np.zeros((L,n))
# Remplissage de Pbis
gamma = np.exp(-r * (tau_star - t)) * f_payoff(S_tau_star)
for i in range(L):
Pbis[i,:] = gamma * Basis[i](S)
Pbis = np.sum(Pbis, axis=1)
assert(np.size(Pbis) == L)
# Compute H
H = Compute_H(S,Basis)
# Inversion de H
H_inv = np.linalg.inv(H)
# Compute alpha
alpha = np.dot(H_inv,Pbis)
return alpha
def Compute_EnsA(Sk,f_payoff,Basis,alpha):
"""
Fonction vérifiant l'inégalité :
f_payoff(Sk) \geq \sum_{l=1}^{L} \alpha_{l} * P_{l}(Sk)
--------
Arguments:
Sk : Valeur de S au temps k donné.
f_payoff : Payoff de l'option étudiée.
Basis : Base fonctionnelle de l'algorithme de Longstaff-Schwarz.
alpha : Coefficients alpha au temps k.
Returns:
A : Booléen (1 si l'inégalité est vraie, 0 sinon)
"""
L = len(Basis)
Somme = 0
for l in range(L):
Somme += alpha[l] * Basis[l](Sk)
A = max(f_payoff(Sk)-Somme,0)
if(A>0):
A = True
return A
def Longstaff_Schwartz_Algo(MatS, T, Basis, f_payoff, r):
"""
Fonction executant l'algorithme de Longstaff-Schwarz à partir de MC réalisations de (S_{t_0},...,S_{t_n}).
--------
Arguments:
MatS : Matrice contenant les MC réalisations de (S_{t_0},...,S_{t_n}). Taille : (MC,n+1)
T : Temps final (Maturité).
Basis : Base fonctionnelle de l'algorithme de Longstaff-Schwarz.
f_payoff : Payoff de l'option étudiée.
r : Rendement de l'actif sans-risque.
K : Strike.
Returns:
B : Prix de l'option bermudienne au temps t0 = 0.
"""
# Initialisation
MC = np.shape(MatS)[0]
n = np.shape(MatS)[1]
B = 0
B_date = np.zeros(n)
# Discrétisation temporelle
t_discret = np.linspace(0,T,n)
# Boucle Monte-Carlo et construction de tau
for j in range(MC):
if(j%100==0):
print(str(j)+"ième itération...")
# Réalisation j
S = MatS[j,:]
tau_star_j = np.zeros(n)
Stau_star_j = np.zeros(np.size(S))
# Initialisation tau
tau_star_j[-1] = T
assert(np.size(t_discret) == n)
# Initialisation alpha
t_interet = t_discret[n-2]
tau_star = tau_star_j[-1]
Stau_star_j[-1] = S[-1]
alpha = Compute_alpha(S,Stau_star_j[-1],Basis,f_payoff,t_interet,tau_star,r)
# Récurrence tau
for p in range(1,n):
# Compute A
A = Compute_EnsA(S[n-1-p],f_payoff,Basis,alpha)
# tau_star maj
tau_star_j[n - 1 - p] = t_discret[n - 1 - p] * A + tau_star_j[n - 1 - p + 1] * (1-A)
# alpha maj
if (p < n-1):
t_interet = t_discret[n - 1 - p - 1]
tau_star = tau_star_j[n - 1 - p]
ind = np.argwhere(t_discret == tau_star)[0][0]
Stau_star = S[ind]
Stau_star_j[n-1-p] = Stau_star #maj Stau_star
alpha = Compute_alpha(S,Stau_star,Basis,f_payoff,t_interet,tau_star,r)
# Prix au temps t = 0
ind = np.argwhere(t_discret == max(tau_star_j))[0][0]
B += np.exp(-r * t_discret[ind]) * f_payoff(S[ind])
B_date += np.exp(-r * tau_star_j) * f_payoff(Stau_star_j)
# Monte-Carlo
B = B/MC
B_date = B_date/MC
return B, B_date
Execution de l'algorithme de Longstaff-Schwarz
# S
t,W = dMouvementBrownien(T,n,MC)
MatS = prixActif(t,W,sg,r,S0)
# Basis
L0 = lambda x : 1
L1 = lambda x : 1 - x
L2 = lambda x : (1/2) * (2 - 4 * x + x**2)
L3 = lambda x : (1/6) * (6 - 18 * x + 9 * x**2 - x**3) #x^3 = x^2*x
Basis = np.array([L0, L1, L2, L3])
B, B_date = Longstaff_Schwartz_Algo(MatS, T, Basis, f_payoff, r)
print("==============")
print("Bermudan Put Price",B)
0ième itération... 100ième itération... 200ième itération... 300ième itération... 400ième itération... 500ième itération... 600ième itération... 700ième itération... 800ième itération... 900ième itération... ============== Bermudan Put Price 9.524233698573683
$\mathbf{Q2})$
fig = plt.figure(figsize=(16,5))
ax = fig.add_subplot();
ax.plot(t,B_date, label="Prix d'exercice Put Bermudien", color='royalblue')
ax.set_title('Comparaison des dates d\'exercice | n = {}, MC = {}'.format(n_input,MC))
ax.set_xlabel('$ \\tau $')
ax.set_ylabel('Prix d\'exercice')
ax.legend()
ax.grid()
plt.show()
fig = plt.figure(figsize=(16,5))
ax = fig.add_subplot();
ax.plot(t[1:],B_date[1:], label="Prix d'exercice Put Bermudien", color='royalblue')
ax.set_title('Comparaison des dates d\'exercice | n = {}, MC = {}'.format(n_input,MC))
ax.set_xlabel('$ \\tau $')
ax.set_ylabel('Prix d\'exercice')
ax.legend()
ax.grid()
plt.show()