Dans ce rapport portant sur les modèles à pics de prix, nous allons tout d’abord présenter chacun des trois différents articles qui composent ce sujet en présentant les modèles et processus qui s’y trouvent ainsi que les utilisations qui en sont faites, et également en discutant de l’estimation des paramètres de ces derniers. Dans un second temps, nous critiquerons ces modèles au vu du contenu du cours, et enfin, dans la dernière partie, nous simulerons des prix sur un des modèles présents dans ces articles, avec une période de simulation libre.
Le but de cet article, réalisé par M.Davison et K.Anderson entre autres, est de construire un modèle stochastique pour les prix spot (établis le jour J pour le lendemain) de l’électricité qui, à la différence d’autres modèles existant, modélise directement le prix, et de l’appliquer au problème de tarification des options sur l’énergie électrique.
La plupart des matières premières sont échangées sur les marchés depuis des années, mais pas l’électricité, d’une part car elle n’est pas stockable, et d’autre part car la plupart des consommateurs sont indifférents au prix à court terme de cette ressource. Cela entraîne donc une grande volatilité des prix, d’où l’existence des pics de prix sur le marché, que l’on observe notamment en été et qui sont l’objet principal de cet article.
Le modèle est une combinaison de deux approches : la première est appelée « top-down ». Elle consiste en l’examination des séries de prix pour estimer les paramètres, ce qui rend le modèle dépendant du marché. La deuxième approche, appelée « bottom-up », comporte des modèles pour la demande et la capacité de génération. On aboutit donc à un modèle hybride, qui permettra de fixer avec précision les prix des contrats à terme et des options sur l’électricité. Il possède 4 facteurs principaux : l’existence des pics de prix, les prix sont souvent nuls hors des pointes, les prix peuvent être parfois négatifs, et les prix de l’électricité ne dérivent pas indéfiniment. Une solution pourrait être de construire un modèle aléatoire dans lequel la série temporelle de prix \(S_t\) serait tirée à partir d’une famille de fonctions de densité de probabilité \(P_t\), une distribution dont les propriétés pourraient varier dans le temps. On souhaite également que cette fonction de densité de probabilité décrivant les prix de l’électricité ait plusieurs pics. Pour cela, on introduit un modèle de commutation (switching model en anglais), qui est un mélange de deux distributions différentes que l’on va décrire par la suite. On a donc la fonction de densité de probabilité des prix sous-jacents suivante: \[P(X)=(1-\epsilon(\alpha(t)))P_L(S)+\epsilon(\alpha(t))P_H(S)\] qui représente le mélange des deux distributions.
Un tel modèle est contrôlé par une variable commutative \(\epsilon\) qui, à chaque pas de temps, détermine de quelle distribution on tire le prix, c’est-à-dire soit la distribution « prix haut » \(P_H(t)\), soit la distribution « prix bas » \(P_L(t)\). Ensuite, à chaque pas de temps encore une fois, on génère une variable aléatoire suivant la loi uniforme sur [0,1], puis on compare la réalisation observée avec le seuil \(\epsilon\) : si elle est inférieure ou égale au seuil, alors il y a un pic de prix et le prix est tiré de \(P_H(t)\), tandis que si elle est supérieure, alors il n’y a pas de pic et donc le prix est tiré de \(P_L(t)\).
De plus, pour bien modéliser les prix de l’électricité et les éventuels pics, il faut prendre en compte les phénomènes inattendus, comme par exemple l’apparition de pics due à une soudaine augmentation de la demande en-dehors des « normes » saisonnières, ou encore lors de restrictions d’approvisionnement. C’est pourquoi notre variable \(\epsilon\) doit être fonction de la demande et de la capacité disponible. Tous ces différents facteurs sont regroupés dans un paramètre adimensionnel \(\alpha\), qui est égal au quotient de la demande instantanée de puissance sur la capacité de génération.
Par conséquent, on peut établir une relation entre \(\epsilon\) et \(\alpha\): en effet, lorsque la capacité excède largement la demande, alors la probabilité d’obtenir un pic est très faible c’est-à-dire que \(\epsilon\) tend vers 0 lorsque \(\alpha\) tend vers 0. Par contre, si toute la capacité est requise pour faire face à la demande, alors on est certain d’obtenir un pic de prix, \(\epsilon\) tend vers 1 lorsque \(\alpha\) tend vers 1.
On a donc la relation fonctionnelle suivante : \[\epsilon(t)=\frac{1}{2}tanh(20^*(\alpha(t)-0.85))+0.5\] Il existe également une valeur limite pour \(\alpha\), qui est environ 0,85. Lorsqu’il est inférieur, les observations de prix élevés sont rares, tandis que si \(\alpha\) est supérieur à 0,85, alors la probabilité d’avoir un pic croît rapidement.
Maintenant que l’on sait cela, l’implémentation de ce modèle peut avoir lieu : on génère un \(\alpha\), à partir duquel on obtient un \(\epsilon\), d’après la formule précédente, qui est la probabilité d’obtenir un pic. Cette variable \(\epsilon\) est ensuite utilisé comme paramètre de commutation pour savoir si l’on doit tirer le prix de la distribution « haute » ou « basse ». Finalement, on aboutit à une série de prix sur une certaine échelle de temps, avec \(P_H\) et \(P_L\) des distributions normales. Le modèle mis en place nécessite également une modélisation de la demande, de la capacité de génération et des prix horaires. Ces derniers sont calculés sur une plage horaire allant de 7h à 23h, et cela donne donc un profil quotidien.
On peut donc résumer ainsi le modèle et ses principales équations, en s’appuyant sur les données existantes sur le marché présent dans l’article, celui de Pennsylvanie, acquises sur 252 jours ouvrés. On a tout d’abord l’espérance de la demande \[E[L(t)]=15000(sin(2\pi(\frac{day}{252})-2)+cos(\pi(\frac{day}{252})+1.5))+40000\] ainsi que la capacité \[C(t)=\left\{
\begin{array}{ll}
58000, & si \; 85<day<200 \\
50000, & sinon
\end{array}
\right.\] d’où la relation suivante incluant \(\alpha\) : \[\alpha(t)=\frac{L(t)}{C(t)}\] Ensuite, étant donné que l’on a déjà la relation entre \(\epsilon\) et \(\alpha\), on peut déterminer si l’on tire dans \(P_H\) ou dans \(P_L\) en comparant la réalisation de la variable uniforme sur [0,1] avec le seuil \(\epsilon\) calculé, et en sachant que les deux distributions suivent les lois normales respectives : \(P_L=N(25,16)\) et \(P_H=N(370,320)\). Et on répète ce processus pour chaque jour de la période considérée, ce qui nous donne la figure suivante :
Passons maintenant aux résultats obtenus, le but étant finalement de tarifier les options sur l’électricité. Pour cela, il a fallu établir un modèle de processus des prix sous-jacents, décrit précédemment, et c’est donc la série de prix quotidiens qui va nous permettre d’estimer la valeur d’une telle option. En utilisant le pricing horaire, les options quotidiennes sont calculées en moyennant ces prix horaires pour donner le prix quotidien moyen entre 7h et 23h, qui sera donc utile pour donner une valeur de l’option (il faut un grand nombre de réalisations). Il faut donc désormais comparer le modèle journalier avec le modèle horaire pour savoir lequel des deux est le plus efficace. Voici donc le tableau suivant :
Pour les prix forward, peu de différences hormis en été. Pour les options Put, pour chaque option avec strike 35, il y a une valeur pour chaque modèle, et l’on peut remarquer un écart entre les deux modèles. Pour l’option Call présente dans le tableau, l’écart entre les deux valeurs est plus important, car le comportement des deux modèles est très distinct, notamment au niveau de l’estimation des paramètres.
En définitive, pour déterminer le meilleur modèle pour tarifier les options, il faut examiner la structure des contrats. L’usage d’un modèle journalier a pour avantage de bien tarifier les options « all-day », alors qu’un modèle horaire donne une meilleure représentation du comportement de hausse des prix, car très sensible aux valeurs de \(\alpha\); en effet, lorsqu’on utilise le modèle journalier, il y a une perte d’informations sur ce comportement. De plus, le modèle horaire permet une plus grande distinction sur les deux distributions. Il en résulte que ce dernier doit être plus approfondi, afin de l’améliorer et de s’en servir à bon escient.
Le but de cet article, réalisé par F.Benth et J.Altyte-Benth, est de modéliser les prix spot sur les marchés de l’énergie avec des processus d’Ornstein-Uhlenbeck exponentiels non-Gaussiens. Il généralise le mouvement Brownien géométrique et le modèle de retour à la moyenne de Schwartz en introduisant les processus de Lévy. De plus, au lieu de modéliser la dynamique des prix spot la solution d’une EDS avec sauts, on s’intéresse à une modélisation statistique, plus avantageuse que la première citée. Cet article discute également du problème des prix à terme et des options.
Le modèle de Schwartz suppose des fluctuations de prix normalement réparties, par conséquent, il semble inapproprié pour modéliser les grandes variations de prix. C’est pourquoi la plupart des auteurs partent de l’EDS qui définit la dynamique de Schwartz et y ajoutent un terme de saut, et ces sauts de prix sont modélisés par des processus de Lévy, dont nous allons parler tout de suite.
Un tel processus (de Lévy-Kintchine) s’écrit : \[L_t = \nu_t + \sigma B_t + \int_{\mid z \mid<1}z\overline{N}((0,t],dz)+\int_{\mid z \mid\ge1}zN((0,t],dz)\] avec \(B_t\) un mouvement Brownien standard, \(\nu\) une constante, \(\sigma\) une constante positive, et N une mesure aléatoire homogène de Poisson avec compensateur \(dtl(dz)\). On a également le processus stochastique des prix spot \(S_t\) : \[S_t=\Lambda(t)e^{X_t}\] avec \(\Lambda\) la saisonnalité, ainsi que la dynamique du processus non-Gaussien d’Ornstein-Uhlenbeck \(X_t\) : \[dX_t=a(m-X_t)dt+dL_t\] où \(X_0=x\), a est la vitesse de retour à la moyenne \(>0\), m une constante positive indiquant la moyenne à long terme du processus. La solution est de la forme : \[X_t=xe^{-at}+\frac{m}{a}(1-e^{-at})+\int_{0}^{t}e^{-a(t-s)}dL_s\]
Lorsque \(X_t=L_t\), c’est-à-dire quand a=0, alors la dynamique des prix spot devient l’exponentielle d’un processus de Lévy, ce qui généralise le mouvement Brownien géométrique classique. Il est aussi question du processus de Lévy gaussien normal inverse (NIG dans l’article), qui permet de modéliser les grandes fluctuations de prix en utilisant la modélisation des résidus des prix spot logarithmiques. Un processus de Lévy est appelé comme ceci si \(L_1\) est distribué selon la distribution NIG, qui dépend de 4 paramètres et qui fait partie de la famille des distributions hyperboliques : \[f(x;\mu,\alpha,\beta,\delta)=\frac{\alpha\delta}{\pi}exp(\delta\sqrt{\alpha^2-\beta^2}+\beta(x-\mu))\frac{K_1(\alpha\sqrt{\delta^2+(x-\mu)^2})}{\sqrt{\delta^2+(x-\mu)^2}}\] avec \(\mu\) un réel indiquant la localisation de la densité, \(\beta\) un réel qui est le paramètre d’asymétrie, \(\delta\) positif qui est un paramètre d’échelle, \(\alpha\) qui mesure l’épaisseur de la queue , et \(K_1\) la fonction modifiée de Bessel du second type et d’ordre 1.
S’en suit une analyse de données sur deux séries de prix. La première montre les prix spot du Brent oil sur 8 ans et la deuxième les prix du SAP (System Average Price) sur 1 an. Pour trouver la saisonnalité annuelle des prix spot, on utilise le développement de Fourier du premier ordre : \[\Lambda(t)=c_1+c_2cos(\frac{2\pi}{n}(t-c_3))\] Les constantes \(c_1\) (le niveau), \(c_2\) (amplitude), et \(c_3\) (décalage temporel) sont calculées en utilisant la méthode des moindres carrés. On peut remarquer pour le SAP des pics beaucoup plus prononcés que pour le Brent oil. Par la suite, les paramètres de la distribution normale (moyenne et écart-type) et de la distribution NIG sont estimés pour chacune des deux séries de prix, et ce pour les rendements logarithmiques des énergies.
En observant les graphes représentant les densités empiriques et ajustées, on se rend compte que la distribution NIG est plutôt bien corrélée avec les queues lourdes.
Maintenant, on s’intéresse au modèle de retour à la moyenne avec \(a>0\). L’estimation se fait en 2 étapes. Tout d’abord, il faut faire une régression linéaire pour estimer la moyenne à long terme ainsi que la vitesse de retour à la moyenne. Ensuite, il faut regarder les propriétés des résidus en les ajustant aux distributions normales et NIG, elles sont observables dans les tableaux 4 et 5 de l’article. La conclusion est que l’hypothèse de normalité est rejetée au seuil 1%.
On s’intéresse maintenant aux prix des produits dérivés sans arbitrage. Il existe un continuum de prix sans arbitrage, un pour chaque mesure de martingale équivalente, sachant que par définition, une mesure de martingale équivalente Q est une mesure de probabilité équivalente à P et tel que tous les actifs échangeables en continu sont des martingales après actualisation. On considère également les modèles de prix spot avec des processus de Lévy plutôt que le mouvement Brownien.
Dans cet article, on dérive des mesures équivalentes par la transformée d’Esscher. Il faut une condition d’intégrabilité exponentielle sur la mesure de Lévy pour s’assurer de l’existence de moments sur le processus de prix. De plus, l’évaluation des forward et des options requiert au moins le moment d’ordre 1, d’où les deux lemmes utilisant l’intégrabilité de la mesure de Lévy qui s’en suivent, avec les preuves, et où il est également question du prix du risque de marché \(\theta\).
L’article aborde ensuite les contrats forward, qui sont une application du modèle, dont le prix à l’instant t pour une livraison à la date \(T>t\) est donné par : \[F^\theta(t,T)=\Lambda(T)exp(\int_{t}^{T}[\phi(\theta(s)+e^{-a(T-s)})-\phi(\theta(s))]ds)(\frac{S_t}{\Lambda(t)})^{exp(-a(T-t))}\] Si l’on dispose de N contrats forward, avec un temps de livraison \(T_i\), \(i=1,...,N\), alors on a : \[F^\theta(0,T_i)=\Lambda(T_i)exp(\int_{0}^{T_i}[\phi(\theta(s)+e^{-a(T_i-s)})-\phi(\theta(s))]ds)(\frac{S_0}{\Lambda(0)})^{exp(-aT_i)}\] Grâce à cela, on peut estimer le \(\theta\) qui minimise l’écart entre les prix théoriques et observés, sous une certaine norme.
Le prix sans arbitrage d’une option avec payoff \(f(S_T)\) à maturité T est donné par : \[P_t=e^{-r(T-t)}E_Q[f(S_T\mid \mathcal{F}_T)]\] Pour trouver le prix en \(t=0\), on doit donc calculer : \[P_0=e^{-rT}E_\theta[f(\Lambda(T)e^{X_T)}]\]
En définitive, cet article a permis d’établir un processus d’Ornstein-Uhlenbeck exponentiel non-Gaussien pour les prix spot des énergies, qui s’ajuste simplement aux données de prix si l’on compare à d’autres modèles stochastiques. De plus, le fait d’utiliser la distribution NIG pour modéliser les résidus dans ce processus permet de rejeter l’hypothèse de normalité et de dire que cette distribution s’ajuste relativement bien aux données et donne une vision plus précise de la dynamique des prix spot, et donc des éventuels pics qui peuvent apparaître.
Cet article, réalisé par A.Cartea et M.Figueroa, présente un modèle de diffusion de saut avec retour à la moyenne pour les prix spot de l’électricité et dérive le forward correspondant, en se basant sur des données spot et forward provenant d’Angleterre et du Pays de Galles.
Comme précisé dans les articles précédents, l’électricité est trop compliquée ou trop chère à stocker. De plus, même si le marché des énergies possède plusieurs similarités avec d’autres marchés, il présente tout de même deux caractéristiques qui le distinguent de ceux-ci : le retour à la moyenne des prix sport et également l’existence de pics de prix. Cela se voit plus particulièrement dans le marché de l’électricité. Ces caractéristiques sont importantes pour ensuite tarifier les produits dérivés de l’énergie. Pour les observer, il existe deux classes de modèle : les modèles spot, et les modèles forward. Cependant, dans certains articles, les modèles spot font apparaître le caractère de retour à la moyenne des prix de l’électricité, mais pas les pics de prix. C’est pourquoi cet article se penche sur cet aspect, en présentant un modèle montrant explicitement ces caractéristiques, ainsi que la saisonnalité. Il s’intéresse aussi à la calibration des paramètres aux marchés anglais et gallois, puis dans un deuxième temps à l’estimation des paramètres du modèle en utilisant à la fois les données spot existantes ainsi que les prix forward actuels.
Passons maintenant à la première partie concernant l’analyse de données sur l’indice FTSE100 et sur les prix moyens quotidiens de l’électricité en Angleterre et au Pays de Galles entre 2001 et 2004.
La figure 1 montre que la trajectoire du prix peut être vue comme la combinaison d’une tendance déterministe et de chocs aléatoires, tandis que sur la figure 2, il y a une forte tendance de retour vers la moyenne, avec de temps à autre des pics. Ces prix de l’électricité ne sont pas pour autant normalement distribués : en effet, la figure 3 nous montre qu’ils ne le sont pas, car on n’observe pas une droite mais plutôt une courbe.
De plus, le modèle de Black-Scholes suppose également une distribution indépendante des rendements, ce qui se voit avec un test d’autocorrélation (le coefficient doit être proche de 0). Cependant, sur les marchés de l’électricité, elle est élevée, ce qui souligne une saisonnalité des rendements (tous les 7 jours d’après la figure 4).
Mais pour estimer les paramètres du modèle, il faut se défaire de celle-ci, en utilisant le rendement désaisonnalisé \(R_t\) défini par : \[R_t=r_t-\overline{r}_d\] avec le rendement au temps t, et la moyenne correspondant au jour que \(r_t\) représente, et on se rend compte sur la figure 5, qui montre le test d’autocorrélation des rendements désaisonnalisés, que la saisonnalité est beaucoup moins perceptible.
En ce qui concerne les pics, un algorithme permettant de retirer les rendements des sauts et donc d’améliorer le test de normalité, sera utile pour calibrer le modèle, en estimant les fréquences cumulées des sauts ainsi que d’autres informations statistiques pertinentes.
Nous pouvons maintenant aborder le modèle présenté par l’article, qui est un modèle basé sur le spot, car il y a un manque de données sur les prix de l’électricité dans les deux pays, ce qui empêche d’utiliser un modèle basé forward. Et, puisque l’on se base sur l’existence du retour à la moyenne et des pics de prix, on peut utiliser un modèle de diffusion de saut à un facteur avec retour à la moyenne, ajusté aux effets de saisonnalité.
Tout d’abord, soit le processus logarithmique des prix suivants : \[lnS_t=g(t)+Y_t\] tel que \(S_t\) peut s’écrire : \[S_t=G(t)e^{Y_t} \; (3)\] avec une fonction déterministe de la saisonnalité et \(Y_t\) un processus stochastique dont la dynamique est : \[dY_t=- \alpha Y_tdt + \sigma(t)dZ_t+lnJdq_t \; (4)\] \(Y_t\), dans (4), est un processus de diffusion de saut à retour moyen de niveau zéro pour le prix sous-jacent \(S_t\), \(\alpha\) la vitesse de retour à la moyenne, \(\sigma(t)\) la volatilité dépendante du temps, J la taille du saut aléatoire proportionnelle, \(dZ_t\) l’incrément du mouvement Brownien, et \(dq_t\) un processus de Poisson tel que : \[dq_t=\left\{
\begin{array}{ll}
1, & avec \; proba \; ld_t \\
0, & avec \; proba \; 1-ld_t
\end{array}
\right.\] avec l la fréquence du processus. De plus, on suppose que J est log-normal. D’après (3) et (4), on peut donc écrire l’EDS vérifiée par \(S_t\) : \[dS_t=\alpha(\rho(t)-lnS_t)S_tdt+\sigma(t)S_tdZ_t+S_t(J-1)dq_t\] où le niveau de retour moyen dépendant du temps est donné par : \[\rho(t)=\frac{1}{\alpha}(\frac{dg(t)}{dt}+\frac{1}{2}\sigma^2(t))+g(t)\] La plupart du temps, \(dqt = 0\), donc on a uniquement le processus de retour à la moyenne. Cependant, à des moments aléatoires, il est possible qu’il y ait un saut, ce qui se traduit par le terme \(S_t-(J-1)\) qui donne le changement avant et après le saut. Cela ressemble au modèle décrit dans le deuxième article, avec des notations différentes.
Il est ensuite question du prix forward. Le prix en t d’un tel produit expirant en T est obtenu comme la valeur attendue du prix spot à l’échéance, sous une mesure de martingale Q équivalente, conditionnellement à l’information disponible en t : \[F(t,T)=E_t^Q[S_t \mid \mathcal{F}_t]\] La suite de cette partie est donc consacrée au calcul de cette espérance, avec une hypothèse forte sur la mesure utilisée : en effet, vu que l’on est sur un marché incomplet, la mesure Q n’est pas unique, mais comme il est compliqué de définir une mesure appropriée, on décide alors de dire que l’on est déjà sous une mesure équivalente. Il suffira par la suite de calibrer le modèle avec des paramètres implicites provenant d’un marché « liquide ». Toutefois, cela reste compliqué pour des marchés récents tels que les marchés de l’électricité anglais et gallois, c’est pourquoi on introduit un prix du risque de marché dans le drift (après avoir appliqué la formule d’Itô à \(ln(S_t))\). Finalement, on arrive à l’expression de prix suivante : \[F(t,T)=G(T) \bigg(\frac{S(t)}{G(t)}\bigg) ^{e^{-\alpha(T-t)}}e^{\int_{t}^{T}\frac{1}{2}\sigma^2(s)e^{-2\alpha(T-s)}-\lambda\sigma(s)e^{-\alpha(T-s)}ds+\int_{t}^{T}\xi(\sigma_J,\alpha,T,s)lds-l(T-t)}\]
où \(\xi(\sigma_J,\alpha\,T,s) \equiv e^{-\frac{\sigma_j^2}{2}e^{-\alpha(T-s)}+\frac{\sigma_j^2}{2}e^{-2\alpha(T-s)}}\)
Nous devons désormais estimer les paramètres du modèle. Une approche est d’utiliser le maximum de vraisemblance, mais elle ne marche pour les données que l’on considère. C’est pourquoi on utilise une approche « hybride » consistant à utiliser les données spot pour estimer la fonction desaisonnalité, la volatilité historique, le taux de réversion moyenne et la fréquence et l’écart-type des pics, et les données de marché forward (à terme) pour estimer le prix du risque de marché.
Tout d’abord, la fonction de saisonnalité est déterministe et elle dépend du marché considéré. Elle peut être sinusoïdale quand un modèle perceptible apparaît entre les mois d’hiver et d’été, ce qui n’est pas évident en Angleterre et au Pays de Galles, ou constante par morceaux dans d’autres cas. Dans notre situation, la fonction correspond aux moyennes mensuelles des données historiques avec une série de Fourier d’ordre 5.
La volatilité, quant à elle, est estimée comme volatilité historique continue, en faisant une moyenne annuelle avec une fenêtre de 30 jours. Le taux de retour à la moyenne est estimé en exécutant une régression linéaire sur les incréments des rendements, et les paramètres de saut sont estimés en utilisant l’algorithme dont nous avons parlé précédemment. Reste le prix du risque de marché, qui est estimé en minimisant les distances au carré des différentes maturités de la courbe forward théorique aux prix du marché d’échéances égales.
Les résultats de ces estimations sont présentés dans la partie suivante.
Enfin, on applique le modèle au prix d’une option d’achat (call) européenne sur un contrat à terme pour le valider. Différentes courbes sont obtenues, comme celles de la simulation de la trajectoire des prix sur un an (figure 8), qui montre la tendance de retour à la moyenne ainsi que l’existence de pics de prix, ou encore les figures 10, 11, et 12 qui représentent des surfaces de prix pour des périodes de temps différentes, respectivement 5 mois, 3 trimestres et 4 saisons, qui montrent également les différents phénomènes que le modèle veut exhiber.
En conclusion, ce modèle a permis d’établir une courbe des prix forward à défaut d’être testé, à cause du manque de données disponibles. Il a aussi permis d’estimer les différents paramètres du modèle (moins nombreux que dans d’autres), grâce aux données spot et forward disponibles, et il permettra de les estimer de façon plus robuste lorsque le nombre de données augmentera. De plus, le résultat obtenu sur la simulation de la trajectoire des prix ressemble à la réalité. Toutefois, des alternatives doivent être étudiées, comme par exemple celles impliquant des processus de Lévy, vus dans le deuxième article.
Les articles que nous avons présentés ci-dessus montre comment les auteurs ont choisi de modéliser les prix de l’électricité spot et forward, en mettant notamment en évidence deux phénomènes importants : le retour à la moyenne et les pics de prix. Nous avons pu voir une modélisation hybride (article 1), qui ne s’appuie pas sur des modèles préexistants, mais plutôt sur un modèle probabiliste avec des distributions de prix aléatoires et des paramètres qui le sont également. Les deuxième et troisième articles étaient quant à eux tournés vers une modélisation à l’aide de modèles et processus connus tels que Lévy ou Schwartz, beaucoup de calculs stochastiques ainsi que l’estimation des paramètres pour calibrer le modèle. Ces deux articles ont également voulu appliquer leur théorie sur des options pour les tarifier et pour montrer la validité de leur modèle, même si chacun estime que des améliorations sont à apporter pour le rendre encore plus performant. En ce qui concerne le cours, certaines approches sont similaires, comme par exemple l’évocation du modèle de Schwartz et les paramètres qui vont avec (saisonnalité, vitesse de retour à la moyenne,…). De plus, il parle de méthodes non citées dans les articles comme celles de Schwartz-Gibson ou Schwartz-Smith, et aussi l’introduction d’autres paramètres tels que le « convenience yield » par exemple, ce qui donne plus de précision à la modélisation des prix selon la ressource choisie. Cependant, comme dit dans la conclusion du cours (slide 123), un des problèmes est l’absence de pics de prix, que les articles tentent de prendre en compte dans leur modélisation. Dans les articles, les auteurs ont fait des choix de modélisation, et tous peuvent être améliorés, notamment lorsque plus de données seront disponibles pour estimer les paramètres, qui deviendront par conséquent plus précis. Ces derniers peuvent toutefois paraître aléatoires, car certains dépendent de phénomènes extérieurs, comme la météo par exemple ou la demande des clients (voir article 1 notamment). Voilà ce que l’on pouvait dire pour ces trois articles, nous pouvons maintenant passer à la simulation des prix, à l’aide des processus de Lévy (fait sur le logiciel R).
Un processus de Lévy est un processus à temps continu, continu à droite limité à gauche, partant de 0, dont les accroissements sont stationnaires et indépendants.
Les exemples les plus connus sont le processus de Wiener ou mouvement brownien et le processus de Poisson
Ainsi dans un premier temps nous avons souhaité simuler un puis des mouvements browniens pour nous familiariser avec les processus de Lévy dans une de leur forme les plus simple et commune.
### Seedage pour avoir toujours le même résultat même si l'on expérimente
### plusieurs résultats au préalable pour voir à quoi cela ressemble
set.seed(2)
### Simulation d’un mvt Brownien
### discretisation du temps
temps = seq(0,1,length=501)
pas.temps = 1/500
### Simulation des accroissements
B.acc = rnorm(500,sd=sqrt(pas.temps))
### Simulation d’une trajectoire
B.sim = c(0,cumsum(B.acc))
plot(temps,B.sim,type="l",xlab="Temps")
## L’option type="l" relie les points entre eux.
Nous avons maintenant simulé 500 trajectoires indépendantes de mouvements Browniens. De plus on s’assure que la trajectoire moyenne empirique est proche de 0.
n.sim=500 ### Nombre de trajectoires
n.point = 201 ### Points de discretisation
temps = seq(0,1,length=n.point)
pas.temps = 1/(n.point-1)
B.acc = matrix(rnorm((n.point-1)*n.sim,
sd=sqrt(pas.temps)),nrow=n.sim)
B.sim = matrix(NA,ncol=n.point,nrow=n.sim)
for (i in 1:n.sim)
{B.sim[i,] = c(0,cumsum(B.acc[i,]))}
dim(B.sim) ## [1] 500 201
## [1] 500 201
B.mean = apply(B.sim,2,mean)
### apply calcule la moyenne pour chaque colonne
plot(temps,B.sim[1,],xlab="temps",type="l",
ylab="Mouvement Brownien",ylim=c(-2,2))
title("Un échantillon de 10 trajectoires")
for (i in 1:10){lines(temps,B.sim[i+1,])}
lines(temps,B.mean,lwd=2,col="red")
legend("bottomleft",c("Trajectoire moyenne"),col=c("red"),lwd=2,cex=0.75)
plot(temps,B.sim[1,],xlab="temps",type="l",
ylab="Mouvement Brownien",ylim=c(-2,2))
title("Un échantillon de 25 trajectoires")
for (i in 1:25){lines(temps,B.sim[i+1,])}
lines(temps,B.mean,lwd=2,col="red")
legend("bottomleft",c("Trajectoire moyenne"),col=c("red"),lwd=2,cex=0.75)
plot(temps,B.sim[1,],xlab="temps",type="l",
ylab="Mouvement Brownien",ylim=c(-2,2))
title("Un échantillon de 50 trajectoires")
for (i in 1:50){lines(temps,B.sim[i+1,])}
lines(temps,B.mean,lwd=2,col="red")
legend("bottomleft",c("Trajectoire moyenne"),col=c("red"),lwd=2,cex=0.75)
Maintenant on s’intéresse à un processus de Lévy plus complexe, en effet on va chercher à combiner un processus de Poisson avec un Processus de Wiener et un drift.
\[L_t = \mu t + B_t + \sum_{i=1}^{N_t}X_i\]
On peut générer un processus de Wiener avec un drift et des lois exponentielles \(W_i\) qui vont correspondre aux durées entre sauts.
Ainsi on a \(N_t = max \left\{ n \in \mathbb{N} : W_1+...+W_n \leq t \right\}\) et les \(W_i\) sont des lois exponentielles indépendantes de moyenne \(\lambda^{-1}\)
### Avec ce code nous générons les 3 composantes
### qui vont nous permettre de créer notre processus
n=1000
h=1/n
lambda=5
set.seed(2)
B=c(0,cumsum(rnorm(n,sd=sqrt(h))))
W=rexp(100,lambda)
N=sum(cumsum(W)<1)
T=cumsum(W[1:N])
X=-rexp(N)
### On crée la fonction du processus
Lt=function(t){
B[trunc(n*t)+1]+sum(X[T<=t])+lambda*t
}
### Enfin on fait le dessin de notre processus
L=Vectorize(Lt)
u=seq(0,1,length=n+1)
plot(u,L(u),type="l",col="blue")
De la même maniere que pour notre mouvement Brownien, on va simuler plusieurs trajectoire du processus de Lévy, afin d’en observer les combinaisons.
### Ici on va créer une matrice afin de stocker nos réalisations
L.mat=matrix(0,500,1001)
for (i in 1:50) {
n=1000
h=1/n
lambda=5
B=c(0,cumsum(rnorm(n,sd=sqrt(h))))
W=rexp(100,lambda)
N=sum(cumsum(W)<1)
T=cumsum(W[1:N])
X=-rexp(N)
Lt=function(t){
B[trunc(n*t)+1]+sum(X[T<=t])+lambda*t
}
L=Vectorize(Lt)
u=seq(0,1,length=n+1)
#plot(u,L(u),type="l",col="blue")
L.mat[i,]=L(u)
}
### Enfin on surperpose nos courbes
### Afin d'obtenir plusieurs trajectoires
L.mean = apply(L.mat,2,mean)
plot(u,L.mat[1,],xlab="u",type="l",
ylab="Processus de Lévy L(u)",ylim=c(-7.5,7.5))
title("Un échantillon de 10 trajectoires")
for (i in 1:10){lines(u,L.mat[i+1,])}
lines(u,L.mean,lwd=2,col="red")
legend("bottomleft",c("Trajectoire moyenne"),col=c("red"),lwd=2,cex=0.75)
### Enfin on surperpose nos courbes
### Afin d'obtenir plusieurs trajectoires
L.mean = apply(L.mat,2,mean)
plot(u,L.mat[1,],xlab="u",type="l",
ylab="Processus de Lévy L(u)",ylim=c(-7.5,7.5))
title("Un échantillon de 25 trajectoires")
for (i in 1:25){lines(u,L.mat[i+1,])}
lines(u,L.mean,lwd=2,col="red")
legend("bottomleft",c("Trajectoire moyenne"),col=c("red"),lwd=2,cex=0.75)
### Enfin on surperpose nos courbes
### Afin d'obtenir plusieurs trajectoires
L.mean = apply(L.mat,2,mean)
plot(u,L.mat[1,],xlab="u",type="l",
ylab="Processus de Lévy L(u)",ylim=c(-10,10))
title("Un échantillon de 50 trajectoires")
for (i in 1:50){lines(u,L.mat[i+1,])}
lines(u,L.mean,lwd=2,col="red")
legend("bottomleft",c("Trajectoire moyenne"),col=c("red"),lwd=2,cex=0.75)
En vertu de ces simulations de processus de Wiener et de processus de Lévy, nous pouvons obtenir des prix dans les modèles définis par nos processus.