Around the Hull-White short-rate model

Apr 30, 2017 11 min read models, project

Introduction
Le modèle de vasicek et quelques pré-requis
Le modèle de Hull et White
Simulations
Conclusion
Annexes

Introduction

Disclaimer : this is the summary of a project done with 2 friends for Pr Jiao Ying in 2017. The corresponding code is avaliable here. Sorry but the following is in french.

Nous allons traiter ici le modèle d’Hull et White, développé à partir du modèle de Vasicek en 1990. Nous commencerons par quelques prérequis, notamment le modèle de vasicek ainsi que des considérations sur les processus d’Ornstein-Uhlenbeck correspondant, avant de nous attaquer au modèle HW en lui-même. Rappelons que nous nous plaçons dans un modèle de marché parfait (infinie liquidité des actifs, coûts de transaction nuls, actifs infiniment divisibles notamment).

Le modèle de vasicek et quelques pré-requis

Rappelons ici la dynamique du taux court dans le modèle de vasicek:

Sous le modèle de vasicek, le taux court suit la dynamique suivante:

$d r_{t} = a (b - r_{t}) d t + σ d W_{t}$

Avec :

$a$ constante positive, représentant la force de rappel
$b$ constante positive, représentant le taux long terme
$σ$ constante positive, représentant la volatilité

On pourrait tout a fait noter $ϑ = a b$ , les paramètres du modèle devenant alors $ϑ$ , $a$ et $σ$

Sous vasicek, le taux court suit un processus d’Ornstein-Uhlenbeck, et la solution correspondante est :

$r_{t} = r_{0} e^{- a t} + b (1 - e^{- a t}) + σ e^{- a t} \int_{0}^{t} e^{a u} d W_{u}$ $O - U$

A t fixé, ce processus d’Ornstein-Uhlenbeck suit une loi gaussienne de paramètres son espérance et sa variances, un processus d’Ornstein-Uhlenbeck est en effet à la fois Gaussien et Markovien.

Soit $X$ une variable aléatoire suivant une loi normale. Alors $e^{X}$ suit une loi log-normale et donc :

$E (e^{X}) = e^{E (X) + \frac{1}{2} V (X)}$

Le modèle de Hull et White

Développé par J.C.Hull et A.White a partir de 1990, le modèle Hull-White est encore très utilisé sur le marché aujourd’hui. L’objectif premier de ce modèle était de correspondre exactement aux données de marché afin de d’obtenir une courbe des taux dans le modèle correspondant exactement à la courbe empirique des taux "réels" fournie par le marché. L’intuition fut la suivante : ajouter un paramètre (non-stochastique) dépendant du temps dans le drift de la dynamique du taux d’intérêt. La mise en oeuvre de ce modèle nécessite une paramétrisation (le paramètre non-stochastique ajouté) à partir des données de marché.

Dans l’extension proposée par Hull-white du modèle de vasicek, c’est le paramètre $ϑ$ qui devient une fonction déterministe du temps :

Sous le modèle Hull-White, le taux court suit la dynamique suivante :

$d r_{t} = (ϑ_{t} - a r_{t}) d t + σ d W_{t}$

Avec :

$a$ constante positive, représentant la force de rappel
$ϑ_{t}$ une fonction déterministe du temps
$σ$ constante positive, représentant la volatilité

Par intégration (via le lemme d’Ito) de la dynamique du taux court sous Hull-White, on obtient: $r_{t} = r_{s} e^{- a (t - s)} + \int_{s}^{t} e^{- a (t - u)} ϑ_{u} d u + σ \int_{s}^{t} e^{- a (t - u)} d W_{u}$

On a: $\begin{aligned} d (r_{t} e^{a t}) & = r_{t} e^{a t} a d t + e^{a t} d r_{t} \\ = r_{t} e^{a t} a d t + e^{a t} ϑ_{t} d t - r_{t} e^{a t} a d t + e^{a t} σ d W_{t} \\ = e^{a t} (ϑ_{t} d t + σ d W_{t}) \end{aligned}$

Ce qui implique que :

$r_{t} e^{a t} = r_{s} e^{a s} + \int_{s}^{t} e^{a u} ϑ_{u} d u + σ \int_{s}^{t} e^{a u} d W_{u}$

et donc :

$r_{t} = r_{s} e^{- a (t - s)} + \int_{s}^{t} e^{- a (t - u)} ϑ_{u} d u + σ \int_{s}^{t} e^{- a (t - u)} d W_{u}$

Sous le modèle Hull-White, le prix d’un zéro-coupon est donné par : $B (t, T) = e x p (\frac{r_{t} (1 - e^{- a (T - t)})}{a} - \frac{1}{a} \int_{t}^{T} ϑ_{u} (1 - e^{- a (T - u)}) d u + \frac{1}{2} [\int_{t}^{T} \frac{σ^{2}}{a^{2}} (1 - e^{- a (t - u)^{2}} d u)])$

On commence par calculer $\int_{t}^{T} r_{s} d s$ directement à partir de la dynamique de r : $\begin{aligned} r_{T} - r_{t} & = \int_{t}^{T} ϑ_{u} - a r_{u} d u + \int_{t}^{T} σ d W_{u} \\ = \int_{t}^{T} ϑ_{u} d u - \int_{t}^{T} a r_{u} d u + σ \int_{t}^{T} d W_{u} \end{aligned}$ Remarquons que grâce a la propriété précédente,

$\begin{aligned} r_{T} - r_{t} & = r_{t} e^{- a (T - t)} + \int_{t}^{T} e^{- a (T - u)} ϑ_{u} d u + σ \int_{t}^{T} e^{- a (T - u)} d W_{u} - r_{t} \\ = - r_{t} (1 - e^{- a (T - t)}) + \int_{t}^{T} e^{- a (T - u)} ϑ_{u} d u + σ \int_{t}^{T} e^{- a (T - u)} d W_{u} \end{aligned}$

Que l’on injectera dans notre équation : $\begin{aligned} a \int_{t}^{T} r_{u} d u & = \int_{t}^{T} ϑ_{u} d u - (r_{T} - r_{t}) + σ \int_{t}^{T} d W_{u} \\ = \int_{t}^{T} ϑ_{u} d u - (- r_{t} (1 - e^{- a (T - t)}) + \int_{t}^{T} e^{- a (T - u)} ϑ_{u} d u + σ \int_{t}^{T} e^{- a (T - u)} d W_{u}) + σ \int_{t}^{T} d W_{u} \\ = r_{t} (1 - e^{- a (T - t)}) + \int_{t}^{T} ϑ_{u} (1 - e^{- a (T - u)}) d u + σ \int_{t}^{T} (1 - e^{- a (T - u)}) d W_{u} \end{aligned}$ Finalement, on a: $\int_{t}^{T} r_{u} d u = \frac{r_{t} (1 - e^{- a (T - t)})}{a} + \frac{1}{a} \int_{t}^{T} ϑ_{u} (1 - e^{- a (T - u)}) d u + \frac{σ}{a} \int_{t}^{T} (1 - e^{- a (T - u)}) d W_{u}$ Calculons son espérance, en notant que l’intégrale sur le brownien est une martingale (d’espérance nulle) : $\begin{aligned} E [\int_{t}^{T} r_{u} d u | F_{t}] & = \frac{r_{t} (1 - e^{- a (T - t)})}{a} + \frac{1}{a} \int_{t}^{T} ϑ_{u} (1 - e^{- a (T - u)}) d u + \underset{= 0}{\underset{⏟}{E [\frac{σ}{a} \int_{t}^{T} (1 - e^{- a (T - u)}) d W_{u} | F_{t}]}} \\ = \frac{r_{t} (1 - e^{- a (T - t)})}{a} + \frac{1}{a} \int_{t}^{T} ϑ_{u} (1 - e^{- a (T - u)}) d u \end{aligned}$ Ainsi que sa variance :

$V a r [\int_{t}^{T} r_{u} d u | F_{t}] = V a r [\underset{d e t e r m i n i s t e}{\underset{⏟}{\frac{r_{t} (1 - e^{- a (T - t)})}{a}}} + \underset{d e t e r m i n i s t e}{\underset{⏟}{\frac{1}{a} \int_{t}^{T} ϑ_{u} (1 - e^{- a (T - u)}) d u}} + \underset{s t o c h a s t i q u e}{\underset{⏟}{\frac{σ}{a} \int_{t}^{T} (1 - e^{- a (T - u)}) d W_{u}}}]$

Et ainsi : $\begin{aligned} V a r [\int_{t}^{T} r_{u} d u | F_{t}] & = V a r [\frac{σ}{a} \int_{t}^{T} (1 - e^{- a (T - u)}) d W_{u}] \\ = \int_{t}^{T} \frac{σ^{2}}{a^{2}} (1 - e^{- a (t - u)})^{2} d u \end{aligned}$

Comme $\int_{t}^{T} r_{u} d u | F_{t}$ est un processus d’Ornstein-Uhlenbeck, on utilise la propriété O-U, et on obtient : $\begin{aligned} B (t, T) & = E [e^{- \int_{t}^{T} r_{s} d s} | F_{s}] \\ = e x p (- E [\int_{t}^{T} r_{u} d u | F_{t}] + \frac{1}{2} V a r [\int_{t}^{T} r_{u} d u | F_{t}]) \\ = e x p (\frac{r_{t} (1 - e^{- a (T - t)})}{a} - \frac{1}{a} \int_{t}^{T} ϑ_{u} (1 - e^{- a (T - u)}) d u + \frac{1}{2} [\int_{t}^{T} \frac{σ^{2}}{a^{2}} (1 - e^{- a (t - u)^{2}} d u)]) \end{aligned}$

Ainsi, nous avons pu complètement définir le prix de n’importe quel zéro-coupon. Pour des raisons de simplification des calculs, nous prendrons les notation suivante :

$\forall t, T, A (t, T) = \frac{1 - e^{- a (T - t)}}{a}$

$\forall t, T, C (t, T) = - \int_{t}^{T} ϑ_{u} A (u, T) + \frac{σ^{2} A (u, T)^{2}}{2} d u$

Et le prix du zéro coupon s’écrit alors :

$B (t, T) = \exp - r_{t} A (t, T) + C (t, T)$

On voit ici que le modèle de Hull-White à un facteur est bien un modèle affine, étant donné que $\exp C (t, T)$ et $A (t, T)$ sont des fonctions déterministes.

Pour en arriver a notre but (construire la courbe des taux), nous allons maintenant déduire du prix du zéro-coupon le prix des taux moyens continus :

$T a u x d^{'} i n t é r ê t c o n t i n u m o y e n$ Par définition, le taux d’intérêt continu moyen s’écrit : $R (t, T) = \frac{- 1}{T - t} l n (B (t, T))$

Ici, on a : $R (t, T) = \frac{- 1}{T - t} (- r_{t} A (t, T) + C (t, T))$

Simulations

Pour l’ensemble des simulations suivantes, nous avons arbitrairement fixé :

$a = 0.5$ , $σ = 1$ , $r_{0} = 0$ et $ϑ (t) = 1 - \frac{\sin (10 t)}{t + 1}$

Comme présenté précédemment, la mise en oeuvre de ce modèle nécessité une paramétrisation de $ϑ$ à partir des données de marchés. Nous choisirons des fonctions simples pour $ϑ$ , sans pousser plus loin la paramétrisation du modèle avec les données de marché.

On fixe le nombre de discrétisation de notre brownien $n$ .
Soit $W$ notre vecteur contenant les valeurs du mouvement brownien
$W_{0} = 0$ et pour $i$ de $1$ à $n - 1$ : $W_{i + 1} = W_{i} + \frac{1}{n} \cdot η$
Avec $η$ une réalisation d’une $N (0, 1)$ pour chaque itération de $i$ .

Simulation d’un mouvement brownien sur 10 ans avec 1 000 000 itérations

Simulation du taux court

Maintenant que nous avons un brownien standard, nous pouvons nous en servir pour simuler le taux court $r_{t}$ .

On fixe $a$ , $σ$ et $ϑ$ .
On génère un mouvement brownien $W$ .
$t m p_{1}$ reçoit le calcul numérique de $\int_{0}^{t} e^{- a \cdot (t - u)} \cdot ϑ (u) d u$ .
On décompose : $\int_{0}^{t} e^{a u} d W_{u} = e^{a t} W_{t} - \int_{0}^{t} a \cdot e^{a u} \cdot W_{u} d u$ .
$t m p_{2}$ reçoit le calcul numérique de $\int_{0}^{t} a \cdot e^{a u} \cdot W_{u} d u$ (méthode des rectangles par exemple).
$r_{t}$ reçoit finalement $r_{0} \cdot e^{- a t} + t m p_{1} + σ e^{- a t} \cdot (e^{a t} \cdot W_{t} - t m p_{2})$ .

Le graphique suivant représente le taux court en fonction du temps $t$ ainsi que de nos paramètres fixés préalablement. Notons que notre fonction $ϑ$ nous permet d’obtenir un taux court qui croît globalement et qui oscille légèrement pour enfin visiblement converger (car $ϑ$ est bornée).

Simulation d’un taux court sur 10 ans avec 1 000 000 de discrétisations du brownien sousjacent

L’intérêt de la fonction sinus dans $ϑ$ était d’avoir un caractère oscillatoire dans notre taux court, l’idée étant de se rapprocher de la réalité sans aller plus loin dans la paramétrisation.

Simulation du prix du zéro-coupon $B (0, T)$

Si $T = 0$ alors $B (0, T) = 1$ .
On déclare $A (t, T, a) = \frac{1}{a} \cdot (1 - e^{- a \cdot (T - t)})$ Comme $r_{0} = 0$ , cela nous épargne le calcul de $t m p_{1} = r_{0} \cdot A (0, T, a)$ .
$t m p_{2}$ reçoit le calcul numérique de $\int_{0}^{T} ϑ (s) \cdot A (s, T, a) d s$
$t m p_{3}$ reçoit le calcul numérique de $\frac{1}{2} \cdot σ^{2} \cdot \int_{0}^{T} A (s, T, a)^{2} d s$
Enfin on a $B (0, T) = e^{t m p_{1} - t m p_{2} + t m p_{3}}$

Le graphique suivant représente donc le prix d’un zéro-coupon $B (0, t)$ en $0$ en fonction de sa maturité $T$ . Notons qu’ici le taux court n’est pas exploité.

Simulation du prix d’un zéro-coupon de maturité $T$

Il est cohérent d’obtenir une fonction décroissante en T en t=0 fixé pour des taux d’intérêt positifs.

Et le suivant représente le prix d’un zéro-coupon $B (t, T)$

Simulation de la courbe des taux B(temps,Temps) sur 10 ans

Simulation et graphe de la courbe de taux

Maintenant que nous sommes capables de calculer le prix du zéro-coupon, il ne reste plus qu’à l’inclure dans le calcul des valeurs de la courbe des taux pour enfin pouvoir en sortir une surface.

$R (0, T) = - \frac{1}{T} \cdot \log (B (0, T))$

Simulation de la surface des taux en $t = 0$ sur 10 ans

Nous obtenons bien une courbe des taux concave. Néanmoins, il se trouve que notre courbe de taux finit par être décroissante à partir de $3$ ans. Cette caractérisation contre intuitive est très probablement due à la paramétrisation de notre modèle.
Enfin, le graphique suivant représente la surface des taux.

Simulation de la surface des taux R(temps,Temps) sur 10 ans

Conclusion

Naturellement, nous n’avons pas pu aller plus loin dans l’analyse financière de nos résultats car nous sommes limités par la sélection des paramètres. Pour calibrer un tel modèle, on pourra citer "Calibration Methods of Hull-White Model" article de Sebastien Gurrieri, Masaki Nakabayashi et Tony Wong parut en 2011.

Annexes

Ci-dessous le code que nous avons utiliser pour générer les graphiques. Nous avons utiliser une méthode de Monte-Carlo pour calculer les intégrales stochastiques du taux court, et les formules fermées démontrées ci-dessus pour calculer les zéro-coupons et la courbe des taux.

## Simulation du brownien standard
T <- 1
n <- 100000
d <- T/n
w <- numeric(n)
for(i in 1:n) w[i+1] <- w[i] + d * rnorm(1)
## Graphique corespondant : 
plot.ts(w)
abline(h = 0)
#Donc W est un brownien sur [0;n]


# parametres du modele : 
.alpha = 0.5
.sigma = 1
.r0 = 0
vega <- function(t)
{
  return(1- sin(10*t)/(t+1))
}


#Fonction calculans le taux court r_t : 
r_t <- function(t,r0 = .r0,alpha = .alpha ,volatilite =1,periode_final=10,wienner=w)
{
  # condition d'arret
  if(t==0){
    return(r0)
  }
  
  # Integration prerequise
  max_iter <- length(w)
  somme_seq <- seq(0, round(max_iter*t/periode_final,0))
  N <- length(somme_seq)
  somme = 0
  for(i in 1:N)
  {
    somme = somme + exp(alpha * i * t / N) * w[i]
  }
  somme <- somme / N
  integrale_pour_wienner <- exp(alpha*t) * w[round(max_iter*t/periode_final,0)] - alpha * somme
  
  # return the result : 
  membre_2 <- integrate(function(u) return(exp(-alpha*(t-u))*vega(u)), lower = 0, upper = t)$value
  membre_3 <- volatilite * exp(-alpha*t) * integrale_pour_wienner
  return(r0*exp(-alpha*(t))+membre_2+membre_3)
}

#On implemente notre H&W
temps <- seq(0,10, length = 1000)
r_seq <- temps
for(i in 1:(length(temps)))
{
  r_seq[i] <- r_t(temps[i])
  print(c(1,i))
}
plot(temps[1:length(temps)], r_seq, type = "l", ylab = "r_t", xlab = "time")








#Fabriquons une courbe de taux
#On commence par creer la fonction qui nous donne les B(0,theta)
#B(0,theta) indexe par alpha

B <- function(t, T, alpha=.alpha, volatilite=.sigma)
{
  #condition d'arret
  if(T==0){return(1)}
  
  #prerequis
  A <- function(t,T,alpha=0.5){
    return((1-exp(-alpha*(T-t)))/alpha)
  }
  
  #lets go
  membre_1 <- r_t(t,0, alpha, volatilite)*A(t,T,alpha)
  membre_2 <- integrate(function(s) return(vega(s)*A(s,T,alpha)), lower = 0,  upper = T)$value
  membre_3 <- 1/2 * volatilite^2 * integrate(function(s) return(A(s,T,alpha)^2), lower = 0,  upper = T)$value
  return(exp(- membre_1 - membre_2 + membre_3))
}

# Maintenant faisont une courbe des taux zero coupon en t=0 d'echeance entre 0 et 10 : 
temps <- seq(0,10, length = 100) #t pouvant aller jusqu'a 10 ans ici
B_seq <- temps
for(i in 1:(length(temps)))
{
  B_seq[i] <- B(0,temps[i])
  print(c(2,i))
}
plot(temps[1:length(temps)], B_seq, type = "l", xlab = "time",  ylab = "B(0,T)")




# Passons a la courbe des taux : 

R <- function(t,T, alpha=.alpha,volatilite=.sigma)
{
  return((-1/T) * log(B(t,T,alpha,volatilite)))
}

# Et la courbe des taux en 0 est : 
temps <- seq(0,10, length = 100)
R_seq <- temps
for(i in 1:(length(temps)))
{
  R_seq[i] <- R(0,temps[i])
  print(c(3,i))
}
plot(temps[1:length(temps)], R_seq, type = "l", xlab = "time", ylab = "Valeur des taux")




## Maintenant, occupons nous de charter les ZC et la courbe des taux sur les deux parametres temps : 
temps <- seq(0,10, length = 100)
Temps <- seq(0,10,length=100)
mat <- matrix(NA,nrow=length(temps),ncol=length(Temps))
for (i in 1:length(temps)){
  for (j in 1:length(Temps)){
    if((temps[i]+Temps[j]) < 10){
      mat[i,j] <- B(temps[i],Temps[i])
      print(c(i,j))
    }
  }
}
persp(temps,Temps,mat,theta=30,phi=30)


## Ensuite, la surface des taux moyens continus : 
temps <- seq(0,10, length = 100)
Temps <- seq(0,10,length=100)
mat2 <- matrix(NA,nrow=length(temps),ncol=length(Temps))
for (i in 1:length(temps)){
  for (j in 1:length(Temps)){
    if((temps[i]+Temps[j]) < 10){
      mat2[i,j] <- R(temps[i],Temps[i])
      print(c(i,j))
    }
  }
}
persp(temps,Temps,mat,theta=30,phi=30)

References :

Azizieh, Celine. 2014. “Modélisation Des Taux d’intérêt.” Université Libre De Bruxelle.

Brigo, Damiano, and Fabio Mercurio. 2006. Interest Rate Models - Theory and Practice: With Smile, Inflation and Credit. 2nd ed. Springer Finance. Springer.

Filipovic, Damir. 2009. Term-Structure Models: A Graduate Course. 1st ed. Springer Finance. Springer-Verlag Berlin Heidelberg.

Musiela, Marek, and Marek Rutkowski. 2010. Martingale Methods in Financial Modelling. 2nd ed. 2005. Corr. 3rd printing. Stochastic Modelling and Applied Probability. Springer.

models rates Black-sholes R

Oskar Laverny

Maître de Conférence

What would be the dependence structure between quality of code and quantity of coffee ?