{ "cells": [ { "cell_type": "markdown", "id": "92218b30-7326-428e-822c-48919e55674b", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "# Ornstein-Ulenbeck, Black-Scholes, initiation `pytorch`" ] }, { "cell_type": "code", "execution_count": null, "id": "3d87b7eb-c418-4d0a-af6b-bca19abbf1dc", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "import numpy as np\n", "import scipy.stats as sps\n", "import pandas as pd\n", "import matplotlib.pyplot as plt\n", "import seaborn as sns\n", "sns.set_theme()\n", "from numpy.random import default_rng, SeedSequence\n", "\n", "sq = SeedSequence()\n", "rng = default_rng(sq);" ] }, { "cell_type": "markdown", "id": "86baf321-08aa-4a6a-b958-5bb5dd62d432", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "## Ornstein-Uhlenbeck" ] }, { "cell_type": "markdown", "id": "fdf18848-1bee-44c2-842d-8ec24794a16c", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "On considère $(X_t)_{t \\in [0,T]}$ processus à valeurs dans $\\mathbf{R}$ solution de l'EDS \n", "\\begin{equation*}\n", " \\operatorname{d}\\! X_t = \\lambda (\\mu - X_t) \\operatorname{d}\\! t + \\sigma \\operatorname{d}\\! W_t, \\quad X_0 = x \\in \\mathbf{R},\n", "\\end{equation*}\n", "où $\\lambda > 0$, $\\mu, \\sigma \\in \\mathbf{R}$.\n", "\n", "On pose $Y_t = e^{\\lambda t} (X_t - \\mu)$ qui vérifie \n", "\\begin{equation*}\n", "\\begin{aligned}\n", " \\operatorname{d}\\! Y_t &= \\lambda e^{\\lambda t} (X_t - \\mu) \\operatorname{d}\\! t + e^{\\lambda t} \\operatorname{d}\\! X_t, \\\\\n", " &= \\sigma e^{\\lambda t} \\operatorname{d}\\! W_t,\n", "\\end{aligned}\n", "\\end{equation*}\n", "et $Y_0 = x-\\mu$. \n", "\n", "Le processus $(Y_t)_{t \\in [0,T]}$ est une intégrale de Wiener, donc un processus de Markov et un processus Gaussien de fonction moyenne $m_Y(t) = x-\\mu$ et de fonction de covariance \n", "\\begin{equation*}\n", "\\begin{aligned}\n", " K_Y(s, t) &= \\operatorname{cov}(Y_s, Y_t) = \\sigma^2 \\mathbf{E}\\bigg[\\int_0^s e^{\\lambda u} \\operatorname{d}\\! W_u \\int_0^t e^{\\lambda u} \\operatorname{d}\\! W_u \\bigg], \\\\\n", " &= \\frac{\\sigma^2}{2 \\lambda} \\big( e^{2 \\lambda \\min(s,t)} - 1 \\big).\n", "\\end{aligned}\n", "\\end{equation*}\n", "En particulier, \n", "\\begin{equation*}\n", " \\forall 0 \\le s \\le t, \\quad \\mathcal{L}\\big(Y_t| \\sigma(Y_u, u \\le s)\\big) \n", " \\sim \\mathcal{N}\\Bigl( Y_s ; \\frac{\\sigma^2}{2 \\lambda} \\bigl( e^{2 \\lambda t} - e^{2 \\lambda s}\\bigr) \\Bigr) \n", "\\end{equation*}\n", "\n", "On en déduit que $(X_t)_{t \\in [0,T]}$ un aussi un processus Gaussien de fonction moyenne $m_X$ et de fonction de covariance $K_X$ définies par \n", "\\begin{equation*}\n", " m_X(t) = \\mu + (x-\\mu) e^{-\\lambda t} \\quad \\text{et} \\quad \n", " K_X(s, t) = \\frac{\\sigma^2}{2 \\lambda} e^{-\\lambda (s+t)} \\big( e^{2 \\lambda \\min(s,t)} - 1 \\big).\n", "\\end{equation*}\n", "En particulier, \n", "\\begin{equation*}\n", " \\forall 0 \\le s \\le t, \\quad \\mathcal{L}\\big(X_t| \\sigma(X_u, u \\le s)\\big) \n", " \\simeq \\mathcal{N}\\Bigl( \\mu + (X_s-\\mu) e^{-\\lambda(t-s)} ; \\frac{\\sigma^2}{2 \\lambda} \\bigl( 1 - e^{-2\\lambda(t-s)} \\bigr) \\Bigr) \n", "\\end{equation*}" ] }, { "cell_type": "code", "execution_count": null, "id": "9b13d5f4-8d01-440e-b6c2-7a4a08c7edb2", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "lambd = 5\n", "mu = 1.\n", "sigma = 0.3\n", "x0 = 3.0\n", "T = 5" ] }, { "cell_type": "markdown", "id": "269d0d0d-2c74-4b6f-afc9-96df3eed786c", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Au temps $T = 5$ on a $X_T$ qui suit une loi Normale de paramètres:" ] }, { "cell_type": "code", "execution_count": null, "id": "ad14e9bc-0d59-4755-b855-b67503778c2f", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "mean = mu + (x0-mu)*np.exp(-lambd*T)\n", "var = sigma**2/(2*lambd)*(1-np.exp(-2*lambd*T))\n", "print(\"mean:\\t\", mean)\n", "print(\"var:\\t\", var)" ] }, { "cell_type": "markdown", "id": "fa39de35-e92f-44b4-bbfe-2ce22afa668e", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "question" ] }, "source": [ "### Question: Simulation via $Y$\n", "\n", "Ecrire une fonction \n", "```\n", " ou_direct_np(noise, delta, x0, lambd, mu, sigma)\n", "```\n", "qui prend en argument `noise` un `np.array` de shape `(n, M)` où `n` est le nombre de pas de temps et `M` le nombre de trajectoires à simuler, `delta` est $\\delta = T/n$ et `x0`, `lambd`, `mu` et `sigma` sont les paramètres du processus d'Ornstein-Uhlenbeck.\n", "\n", "L'implémentation doit être la suivante: on simule les accroissements de $(Y_t)_{t \\in [0,T]}$ aux instants $t_k = k \\delta$, puis on reconstruit les $Y_{t_k}$ (via `np.cumsum`) et on obtient $X_{t_k}$ à partir de $Y_{t_k}$." ] }, { "cell_type": "code", "execution_count": null, "id": "c5e04807-9c2e-4401-9016-139c262a29dc", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "aremplir" ] }, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "b3d399bf-9cdf-45ae-8e43-6f163d315f36", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "question" ] }, "source": [ "### Question: Simulation récursive\n", "\n", "De façon similaire, écrire une fonction \n", "```\n", " ou_rec_np(noise, delta, x0, lambd, mu, sigma)\n", "```\n", "qui utilise une boucle `for` pour construire itérativement $X_{t_{k+1}}$ à partir de $X_{t_k}$ (les $M$ trajectoires en même temps)." ] }, { "cell_type": "code", "execution_count": null, "id": "603328a7-867f-4929-a557-7feceb08b483", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "aremplir" ] }, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "08722d69-9152-48c0-b28a-28d8e5e07a3a", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "question" ] }, "source": [ "### Question: vérification code et visualisation\n", "\n", "Vérifier en utilisant le même argument `noise` que les deux fonctions renvoient les mêmes trajectoires. Vérifier aussi la moyenne empirique et la variance empirique terminale.\n", "\n", "Tracer 100 trajectoires." ] }, { "cell_type": "code", "execution_count": null, "id": "fc02b2be-a9db-4e02-9042-375d14ae07ad", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "aremplir" ] }, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "47ea1cc5-b248-4f37-af1c-8b4cb7b59ce7", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "question" ] }, "source": [ "### Question: comparaison des temps \n", "\n", "Comparer les temps de calcul pour différentes valeurs de `n` et `M`." ] }, { "cell_type": "code", "execution_count": null, "id": "9fb70cf2-1d67-421c-9338-333344c766df", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "aremplir" ] }, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "221ddcf9-9740-445b-a49e-c1c2100b9c21", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "## Algorithmes en `pytorch`\n", "\n", "Reprendre les deux fonctions précédentes pour adapter le code à `pytorch`: on remplace les `np.array` de `numpy` par des `torch.tensor` de `pytorch`.\n", "\n", "[Documentation officielle à lire.](https://pytorch.org/tutorials/beginner/basics/tensorqs_tutorial.html)" ] }, { "cell_type": "markdown", "id": "a7442522-29d3-4366-b800-aef10e70d733", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "question" ] }, "source": [ "### Question: Implémentation sur CPU\n", "\n", "Ecrire des fonctions `ou_direct_torch` et `ou_rec_torch` pour utiliser `pytorch` à la place de `numpy`.\n", "\n", "Il suffit d'importer le module `torch` et de remplacer les appels de `numpy` par les fonctions/objets équivalents en pytorch. " ] }, { "cell_type": "code", "execution_count": null, "id": "c8916722-d616-42c1-b6ec-471f21a4a464", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "aremplir" ] }, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "5a1f2faa-15d7-4ed2-be7d-fb8b5a2cbf06", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "question" ] }, "source": [ "### Question: Tester et comparer les temps \n", "\n", "Vérifier vos codes en calculant moyenne et variance sur l'échantillon au temps final. \n", "\n", "Mesurer les temps de calculs et les comparer avec `numpy`." ] }, { "cell_type": "code", "execution_count": null, "id": "c4fe6dc8-b6a0-486b-a8ce-872713f4c5af", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "aremplir" ] }, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "8ef7f6f0-2556-4a9e-902c-d4000f2f43d1", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "question" ] }, "source": [ "### Question: Implémentation sur un autre device\n", "\n", "Pour utiliser un autre device, il faut modifier son code pour demander la création des tensors sur le device en question. \n", "\n", "Adapter le code précédent en ajoutant l'option `device` aux fonctions `ou_direct_torch` et `ou_rec_torch`." ] }, { "cell_type": "code", "execution_count": null, "id": "17db14c7-3c37-4586-8950-9f7600374c28", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "aremplir" ] }, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "e9af6fea-8f2c-48c6-af91-00c37986d05f", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "## Brownien géométrique, modèle de Black-Scholes" ] }, { "cell_type": "markdown", "id": "8c8cacf2-70ca-4d1e-938f-474add4e4a40", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "On considère $(S_t)_{t \\in [0,T]}$ solution de l'EDS sur $[0,T]$\n", "\\begin{equation*}\n", " \\operatorname{d}\\!S_t = r S_t \\operatorname{d}\\!t + \\sigma S_t \\operatorname{d}\\!B_t, \\quad S_0 = x\n", "\\end{equation*}\n", "c'est à dire \n", "\\begin{equation*}\n", " S_t = x \\exp \\bigl((r- \\sigma^2/2) t + \\sigma B_t \\bigr).\n", "\\end{equation*}" ] }, { "cell_type": "markdown", "id": "a9005ed5-0ffb-44b9-b90b-8d17d7b68a04", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "question" ] }, "source": [ "### Question: Simulation directe\n", "\n", "Réfléchir à deux algorithmes différents: l'un en utilisant la transformation exponentielle appliquée à un Brownien (dilaté et translaté), l'autre en utilisant un produit cumulé `cumprod`.\n", "\n", "Ecrire un code de simulation `numpy` puis `pytorch` pour simuler `M` trajectoires aux instants $t_k = k T/n$, $n \\ge 1$." ] }, { "cell_type": "code", "execution_count": null, "id": "55ff1a76-51f7-49c0-9d6c-510a3fa9d8ac", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "aremplir" ] }, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.6" } }, "nbformat": 4, "nbformat_minor": 5 }