Copenhagen Business School 2015
Kandidatafhandling Cand.merc.(mat)
Prisfastsættelse af Gaslagre ved Least Squares Monte Carlo Gas Storage Valuation by Least Squares Monte Carlo Forfattere:
Vejledere:
Frederik Ørnsholt Regli Rasmus Lind
Dato for a evering: 27. maj 2015 Antal normalsider: 120
Nina Lange Kristian Miltersen
Forord Vi vil gerne takke Nina Lange og Kristian Miltersen for deres engagement og vejledning gennem hele processen. Der skal også lyde en tak til Nationalbanken, specielt Morten Kjærgaard, for at have stillet to specialepladser til rådighed. Sidst vil vi gerne takke Mats Broden og Richard Palmieri fra DONG Energy for værdifulde inputs til opgaven.
Abstract Our master thesis consists of two parts. In the rst part we set up a spot price model on natural gas. Like in Müller et al. [2015] and Stoll and Wiebauer [2010] we show that oil prices and temperatures are important factors in uencing the spot price. We extend their model by adding a statistically signi cant eect from gas storage levels. Low storage levels ceteris paribus leads to higher spot prices. In the second part we value a gas storage contract using the Least Squares Monte Carlo approach. We test how the valuation depends on the type of basis functions and various combinations of state variables. We assess the risk associated with a trading strategy purely on the spot price and show how starting the trading strategy in an intrinsic static hedge will diminish the risk.
1
Indhold 1 Indledning
5
2 Problemformulering
7
3 Placering af opgaven i forhold til litteraturen
8
3.1
Litteratur om gasprismodeller . . . . . . . . . . . . . . . . . . . . . . . . . .
8
3.2
Litteratur vedrørende prisfastsættelse af et gaslager . . . . . . . . . . . . . .
10
3.3
Videnskabsteoretisk metode . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
I Spotprismodel
14
4 Databehandling
15
4.1
Gaspriser
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
15
4.2
Temperaturdata . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
16
4.3
Oliedata . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
17
4.4
Gaslagerdata . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
18
5 Gasprismodel
20
5.1
Stoll og Wiebauers gasprismodel . . . . . . . . . . . . . . . . . . . . . . . . .
20
5.2
Müllers udvidelse med olie . . . . . . . . . . . . . . . . . . . . . . . . . . . .
21
5.3
Inddragelsen af gaslagerdata . . . . . . . . . . . . . . . . . . . . . . . . . . .
24
5.4
Temperaturchok . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
25
5.5
Opstilling af model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
26
5.6
Modellering af XtG . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
29
6 Simulationsstudie 6.1
6.2
32
Temperaturmodel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
32
6.1.1
Simulation af Temperatur . . . . . . . . . . . . . . . . . . . . . . . .
34
Olie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
35
2
6.2.1 6.3
6.4
Simulation af oliepriser . . . . . . . . . . . . . . . . . . . . . . . . . .
38
Gaslager model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
39
6.3.1
Simulation af gaslagerbeholdning . . . . . . . . . . . . . . . . . . . .
42
Simulation af gaspriser . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
42
7 Backtesting og diskussion af model 7.1
44
Backtesting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
44
7.1.1
Out of sample på lageråret 14/15 . . . . . . . . . . . . . . . . . . . .
47
7.2
Kausalitet . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
49
7.3
Reduktion i antallet af stier grundet ikke stationær oliemodel . . . . . . . . .
52
8 Kalibrering af spotsimulationer til forwardpriserne
53
8.1
Udglatning af forwardkurven . . . . . . . . . . . . . . . . . . . . . . . . . . .
53
8.2
Modi cerede simulationer . . . . . . . . . . . . . . . . . . . . . . . . . . . .
58
II Prisfastsættelse af en gaslagerkontrakt ved Least Squares Monte Carlo 59 9 Prisfastsættelse af gaslagerkontrakter
60
9.1
Gaslagerkontrakt . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
60
9.2
Forsimplende antagelser . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
62
9.3
Diskretisering af handlingsrum, tid og volume . . . . . . . . . . . . . . . . .
63
9.4
Least Squares Monte Carlo . . . . . . . . . . . . . . . . . . . . . . . . . . . .
68
9.4.1
Basisfunktioner . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
70
Illustrativt eksempel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
73
9.5.1
85
9.5
Rentens påvirkning på LSMC . . . . . . . . . . . . . . . . . . . . . .
10 Sammenkobling af spotprismodellen og LSMC
89
11 Resultater
89
11.1 Basisfunktioner . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
90
11.2 Antallet af simulerede pristier . . . . . . . . . . . . . . . . . . . . . . . . . .
91
11.3 Precision . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
92
11.4 Tilstandsvariable . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
93
11.5 Bias i prissætningen
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
96
11.6 Resultater for de realiserede stier . . . . . . . . . . . . . . . . . . . . . . . .
98
11.7 Hedge og risikomål . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 11.7.1 Statisk hedge ved intrinsic-strategien . . . . . . . . . . . . . . . . . . 101 11.7.2 Risikomål . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104 11.7.3 Risikomål benytte på LSMC værdiansættelsen . . . . . . . . . . . . . 105 11.7.4 Lagerårene 2012-2013 og 2013-2014 . . . . . . . . . . . . . . . . . . . 106 11.7.5 Lageråret 2014-2015 . . . . . . . . . . . . . . . . . . . . . . . . . . . 112 11.7.6 Betydningen af reduktionen af stier . . . . . . . . . . . . . . . . . . . 115 11.8 Opsamling af hovedresultater . . . . . . . . . . . . . . . . . . . . . . . . . . 118
12 Konklusion
120
13 Bilag
121
13.1 Prismodel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121 13.2 2014-2015 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122 13.3 Resultater . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122 13.4 Kode til prismodellen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 13.5 Kode til at danne resultater . . . . . . . . . . . . . . . . . . . . . . . . . . . 140
4
1
Indledning
Naturgas benyttes primært i produktionen af strøm og varme. Efterspørgslen efter gas er præget af sæsonudsving med højere efterspørgsel i vintermånederne end i sommermånederne. Det skyldes den øgede varme- og strømproduktion om vinteren. Priser på naturgas bestemmes i en ligevægt mellem udbud og efterspørgsel, hvor et øget udbud vil medføre et prisfald og en øget efterspørgsel alt andet lige vil give en prisstigning. Den øgede efterspørgsel om vinteren skaber således et spænd i priser mellem sommer og vinter. Dette spænd kan man tjene penge på, såfremt man kan lagre gassen fra sommer til vinter, da man vil være i stand til at købe gassen billigt og sælge den dyrt. Gaslagre spiller derfor en tiltagende vigtigere og vigtigere rolle for udbudssidens sæsonudglatning og imødekommelse af efterspørgslen. I perioden den 1. januar 2012 til den 30. april 2015 er den angivne maksimale tekniske lagerkapacitet i de 28 europæiske lande steget fra 56.244 mcm til 92.271 mcm1 . Lagerkapaciteten er steget markant, men spændet mellem priserne i vintermånederne og sommermånederne er løbende blevet formindsket. Dette kan ses i gur 1.1, hvor SY står for Storage Year. Et Storage Year løber fra 1. april til 31. marts.
Figur 1.1: Kilde: GSE publikationen "The Value of Storage" 20. Marts, 2014. Tilgængelig via http://www.gie.eu/index.php/publications/gse Mindskelsen kan i høj grad tilskrives den forøgede lagerkapacitet. Det betyder også, at 1 Kilde:
GIE
5
lagrene ikke har det samme indtjeningspotentiale som i tidligere år, hvor spændet mellem sommer og vinter var højere. Dette er en udfordring for de kommercielle ejere af gaslagre. Den kommercielle ejer kan enten acceptere den lavere indtjening, der er forbundet med at handle i spændet på forwardkurven, eller kan alternativt løbende handle den mere volatile dagspris. Risikopro len for de to tilgange er vidt forskellige. Ved at handle sæsonspændet i forwardkurven kan man risikofrit låse en gevinst fast, mens det at handle fra dag til dag er behæftet med en del risiko. Det skyldes den løbende stillingtagen til, hvorvidt der skal købes eller sælges. Risikoen i at handle dagligt ligger i, at der kan opstå en uventet udvikling i prisen, hvilket fx kan være et fald i vintermånedernes priser, der så betyder, at ejeren af lageret må sælge sin oplagrede gas med tab. Denne problemstilling for den kommercielle ejer af et gaslager, giver anledning til problemformuleringen i følgende afsnit.
6
2
Problemformulering
I denne afhandling adresserer vi følgende hovedspørgsmål.
• Er det pro tabelt at benytte Least Squares Monte Carlo til at opstille en
handelsstrategi på spotprisen? Opgaven består af to dele. For at besvare hovedspørgsmålet fremsætter vi i del 1 en spotprismodel for naturgas. Her vil vi besvare følgende spørgsmål.
• Hvilke faktorer har en indvirkning på gasprisen?
Specielt vil vi se på eekten fra lagerbeholdningen af gas i EU. • Hvorfor har faktorerne indvirkning på gasprisen? • Hvordan påvirker faktorerne gasprisen? I del 2 vil vi implementere Least Squares Monte Carlo og benytte algoritmen til at prisfastsætte gaslagerkontrakter. Prisfastsættelsen foretages på baggrund af simulationer fra spotprismodellen i del 1. Dette leder til følgende spørgsmål.
• Hvor godt klarer en beslutningsregel, som er baseret på Least Squares Monte Carlo, værdiansættelsen af en gaslagerkontrakt, sammenholdt med en risikofri strategi?
• Hvordan kan risikoen karakteriseres for Least Squares Monte Carlo handelsstrategien på spotmarkedet?
I den forbindelse vil vi kigge på risikomålene Value at Risk og Expected Shortfall.
7
3
Placering af opgaven i forhold til litteraturen
3.1 Litteratur om gasprismodeller Dette afsnit vil bestå af en kort gennemgang af de gængse modeller fra litteraturen, som kan bruges til at modellere gasprisudviklingen. Formålet er at give et billede af, hvad der tidligere er beskrevet i litteraturen, og hvordan vores prismodel placerer sig i forhold til det. Den gængse tilgang til modelleringen af råvarepriser er at benytte faktor modeller enten på spotprisen eller på hele forwardkurven. Der udvælges således en eller ere stokastiske faktorer. Faktor modellerne prøver at inkorporere elementer fra følgende traditionelle teorier:
• The theory of normal backwardation fremsat af Keynes [1930]. Teorien er her, at risikoafdækkende agenter primært optræder på salgssiden i forwardmarkedet. De ønsker således at fastlåse den pris, som de kan sælge til. De er derfor villige til at tilbyde spekulanter et nedslag i salgsprisen. Dette betegnes, som en risikopræmie, da det er en kompensation til spekulanten for at påtage sig risikoen. Derved er forwardprisen lavere end spotprisen, hvilket betegnelsen backwardation dækker over. Det viser sig dog, at forwardkurven ikke altid er i backwardation, så hovedpointen er, at risikoafdækkende agenter tilbyder spekulanter en risikopræmie.
• The theory of storage and convenience yield som af afdækket af Kaldor [1939], Kaldor [1940], Working [1949], Telser [1958] og Brennan [1958]. Her er teorien, at den arbitrage-frie futures pris skal være lig omkostningerne forbundet med at købe og lagre den underliggende råvare minus den convenience yield, der er forbundet med at have råvaren på lager. Convenience yield skal forstås som den nytteværdi, der er forbundet med at have godet. Nytteværdien består i, at man ved at have godet kan forbruge det, hvilket man ikke kan, hvis man først modtager godet på et fremtidigt tidspunkt.
• The Samuelson eect er betegnelsen for det fænomen, at volatiliteten er mindre for forward kontrakter, desto længere tid de har til udløb (Samuelson [1965]).
8
De faktorer, som fortrinsvis benyttes, er spotprisen, convenience yield, renten og langsigtsniveauet for prisen. Den simpleste form er en faktor modeller. Eksempler på en faktor modeller er fx modellen fra Brennan and Schwartz [1985], hvor spotprisen følger en geometrisk Brownsk bevægelse eller Cortazar and Schwartz [2003], hvor en mean reverting process benyttes. I Schwartz [1997] inddrages convenience yield og renten som ekstra faktorer. Schwartz and Smith [2000] opdeler spotprisen i en kort- og langsigtsfaktor, hvilket Boogert and de Jong [2011] udvider med en sommer-vinter faktor til at tage højde for sæsonudsvinget i gassens forwardkurve. Ikke alle disse modeller er fremsat med gasmarkedet for øje. De er i højere grad udtryk for nogle overordnede betragtninger om modelleringen af råvarepriser. En fællesnævner for de ovenstående modeller er, at de bliver estimeret udelukkende på baggrund af informationen, der er indeholdt i forwardkurven. Forwardkurven må formodes at afspejle forventningen til prisens sæsonudsving og vil derigennem latent indeholde information omkring relevante risikofaktorer. Det må dog forventes, at forwardkurven ikke i samme grad er i stand til at beskrive de kortsigtede uktuationer i gasprisen, der opstår ved udbudseller efterspørgselsstød. Information, der relaterer sig til disse stød, må i højere grad være indeholdt i tidsrækkerne for temperatur og beholdningen af gas i de europæiske lagre. Dette skyldes dels, at gasprisen, som nævnt i indledningen, primært bruges i produktionen af varme og strøm, samt at lagerbeholdningen giver en indikation af sikkerhedsnettet på udbudssiden. Sikkerhedsnettet er vigtigt, hvis forbruget skulle stige som følge af enten en lang og kold vinter eller et fald i udbuddet som følge af et stop i eksporten af gas. Dette er tidligere set ved en række geopolitiske kriser. Prismodellen, som opstilles i dette speciale, er inspireret af fremgangen i to artikler. Den første artikel er Stoll and Wiebauer [2010], hvor temperatur anvendes som en eksogen faktor. Mere præcist benyttes normaliserede heating degree days til at forklare gasprisen. Den anden artikel er Müller et al. [2015], hvor der tages udgangspunkt i modellen fra Stoll and Wiebauer [2010]. Olieprisen tilføjes som en eksogen stokastisk faktor, der påvirker gasprisen igennem et 180 dages moving average. Denne tilføjelse skyldes, at gaskontrakter i udpræget grad er indekseret til olieprisen. Med-
9
tagelsen af olieprisen forbedrer ttet markant. Det nye i denne afhandling er at inkludere lagerniveauet i de 28 europæiske lande som en faktor. Idéen er her, at lagerniveauet angiver, hvor velpolstret udbudssiden er og derved angiver evnen til at imødekomme perioder med enten en forhøjet efterspørgsel eller et formindsket udbud. Informationen om lagerdata er oentligt tilgængeligt via Gas Infrastructure Europe, hvor indrapporteringen er gået fra at være ugentlig til daglig. Idéen om at benytte lagerinformation er ikke ny. Dataen fra Gas Infrastructure Europe er således benyttet i Kremser and Rammerstorfer [2012], hvor målet er at modellere den implicerede convenience yield. Til dette bruges lagerinformationen på et ugentligt niveau. Volmer [2011] bruger ligeledes lagerinformation i modelleringen af convenience yield, hvor der benyttes månedligt data fra International Energy Agency. Modellen i denne afhandling adskiller sig fra tilgangene i Volmer [2011] og Kremser and Rammerstorfer [2012] på mindst to punkter. For det første er convenience yield ikke indeholdt i vores model. Vi kan således ikke udtale sig om, hvad convenience yield er. For det andet er modellen en daglig model i diskret tid, hvilket er vigtigt for den senere Least Squares Monte Carlo prisning af lageret. Det er vigtigt, da Cartea et al. [2015] redegør for, at værdien af en lagerkontrakt stiger med antallet af mulige handelstidspunkter. Det er således til formålet en fordel at have en daglig fremfor en ugentlig eller månedlig model.
3.2 Litteratur vedrørende prisfastsættelse af et gaslager Dette afsnit skal danne et kort overblik over litteraturen omkring prissætningen af gaslagre. Værdien af lageret hænger uløseligt sammen med den handelsstrategi, der anvendes for lageret. I litteraturen skelnes der i store træk mellem tre handelsstrategier. Den korte oversigt, som gives her, gengiver i hovedtræk oversigten i Breslin et al. [2008]. De tre metoder er:
• Intrinsic Value er værdien, der kan opnås ved at optimere handlen i spændene for forwardkurven mellem sommer og vinter. Optimeringen sker under bibetingelse af, at der ikke brydes med lagerets tekniske begrænsninger. Denne værdi kan beregnes uden at modellere udviklingen i spotprisen og forwardkurven. Dette baserer sig på, at værdien ikke længere er eksponeret overfor udviklingen i spot- eller forwardprisen, 10
når positionerne i forwardkontrakterne er valgt. Derfor vil vi i denne opgave benytte betegnelsen,
at man indgår et statisk hedge på lige fod med, at man følger intrinsic
strategien. • Rolling Intrinsic Value er en udvidelse af handelsstrategien for intrinsic value. Der foretages en løbende vurdering af, hvorvidt ændringerne i forwardkurven gør det pro tabelt at lukke nogle af de nuværende positioner og påtage sig nye positioner. Værdiansættelsen af lageret for et helt år kommer her til bestå af intrinsic value plus værdien af den løbende rebalancering. Hvad værdien af rebalanceringen bliver vides ikke med sikkerhed på dagen for værdiansættelsen. Denne skal således estimeres, hvilket fx kan ske ved Monte Carlo simulering. Pointen er her, at værdiansættelsen nu er behæftet med usikkerhed forstået på den måde, at den realiserede værdi af lageret ikke nødvendigvis er den samme som, hvad lagerets pris er fastsat til, hvilket var tilfældet for intrinsic strategien. For en mere detaljeret gennemgang af rolling intrinsic strategien henvises til Breslin et al. [2009].
• Spot Optimization baserer sig, som navnet angiver, på at optimere den daglige handel på spotmarkedet. Dette sker typisk ved at opstille en model for spotprisen og så enten bruge en nite dierence tilgang med trinomialtræer eller en simulationsbaseret tilgang ved Least Squares Monte Carlo. I denne afhandling vil den sidstnævnte være i fokus. LSMC blev først fremsat i Longsta and Schwartz [2001], hvor metoden blev brugt til prisfastsættelse af amerikanske optioner. Første gang metoden blev anvendt til prisfastsættelsen af gaslagerkontrakter var i Boogert and De Jong [2008], hvor metoden blev tilpasset på en sådan måde, at den tager højde for lagerets fysiske begrænsninger. Udelukkende at handle på spotmarkedet er mere risikabelt end de to intrinsic strategier, da man køber med en forventning om, at man kan sælge dyrere i fremtiden. Omstændighederne kan dog arte sig sådan, at man ikke er i stand til sælge med gevinst.
11
I denne afhandling sammenholdes tre strategier. Den første er intrinsic strategien, den anden den rene LSMC spotstrategi og den sidste er en blanding af de to strategier. LSMC strategierne baserer sig på simulationerne fra den opstillede spotprismodel i afsnit 5. Spotprismodellens simulationer kalibreres på sidste dag inden prisfastsættelsen til forwardkurven, hvilket er beskrevet i afsnit 8. I denne afhandling opstilles ikke en model, som beskriver udviklingen i forwardpriserne. Der ses således ikke på værdiansættelsen ved rolling intrinsic strategien, da en modellering af udviklingen i forwardprisen så netop ville være en nødvendighed.
12
3.3 Videnskabsteoretisk metode Metoden i denne opgave vil generelt basere sig på den positivistiske tankegang. Positivismen bygger på empirismen, det vil sige på induktion, indsamling af kvantitative data og disses ordning gennem logiske principper. Positivisten mener, at alt kan reduceres til matematisk fysik (fysikalisme). Videnskab er sand, hvis den korresponderer med virkelige fænomener. Årsager og sammenhænge veri ceres ud fra empiri og kan spores tilbage hertil. I opgaven opstilles en statistisk model for gasprisen, der baserer sig på teorier om sammenhænge på markedet. Den historiske gaspris forsøges forklaret ved at benytte olieprisen, temperaturen og lagerniveauet. Gennem den opstillede matematiske model vil der ligeledes estimeres forventede fremtidige gaspriser. Kritikken af den positivistiske tilgang kommer fra ere sider. Socialkonstruktivisterne kritiserer blandt andet positivisternes ordning af data. Positivisterne hævder, at deres indsamling af data er objektiv. Her vil socialkonstruktivisterne indvende, at de har udvalgt de faktorer, der måles på. Der er således foretaget en udvælgelse af, hvad der skal observeres, og herved samtidig en udvælgelse af, hvad der ikke skal observeres. Denne pointe er værd at tage med i modelarbejdet. I modelarbejdet har vi valgt at benytte en række målbare faktorer, som påvirker gasprisen, mens andre faktorer er udeladt. Nogle eksempler på faktorer, der ikke er medtaget, er kulprisen og det geopolitiske forhold mellem gaseksportører og -leverandører. Kulprisen har betydning, da den fungerer som substitut for gassen i produktionen af strøm. Man kunne således have medtaget den som input. De politiske forhold eller nansiel uro på gasmarkedet tillægges også at have en eekt. Eekten af disse er svære at kvanti cere. Det er ligeledes svært at forudsige, hvorvidt en geopolitisk krise opstår, og vil de geopolitiske kriser, der tidligere er set, have den samme eekt som fremtidige geopolitiske kriser? Med kritikken af den positivistiske tilgang i mente, vil vi i denne afhandling bestræbe os på at anskueliggøre empiriske sammenhænge. Vi vil samtidig kritisk diskutere de fundne resultater, da de skal ses i lyset af den subjektive udvælgelse af faktorer, vi har foretaget.
13
Del I
Spotprismodel
14
4
Databehandling
I dette afsnit beskrives data for gasprisen, temperaturen, olieprisen og gaslageret. Information omkring datakilder, der benyttes, og valget af disse er opgivet. Det meste af data bliver behandlet, og diverse skridt i databehandlingen vil blive beskrevet nedenfor.
4.1 Gaspriser Data er fra Title Transfer Facility (TTF) i Holland, da dette er en de største trading hubs i Europa. En tidsrække for gasprisen stykkes sammen af to tidsrækker. Den første tidsrække er day ahead prisen, som er prisen på gas med levering dagen derpå. Denne er hentet fra Bloomberg med koden TTFGDAHD OECM. Disse priser er kun opgivet på arbejdsdage. Den anden tidsrække er weekendprisen, der ligeledes er hentet fra Bloomberg med koden TTFGWKND OECM. Over weekenden er der ingen priser til rådighed, så den mest retvisende pris i weekenden må alt end lige være den sidste weekendpris, der er oentliggjort. Denne fremskrives, så prisen på gas lørdag og søndag er den samme. Det samme er gældende for helligdage. Den samlede tidsrække vil i opgaven benævnes som spotprisen. Stoll and Wiebauer [2010] tager højde for outliers. Dette skyldes, at helt ekstraordinære begivenheder kan påvirke gasprisen. En gaskrise mellem Ukraine og Rusland, som de har haft ere af i det nye årtusind, kan have stor påvirkning på gasprisen, men gaskriser kan samtidig være svære at modellere. Data, der er korrigeret for outliers, benyttes ikke igennem hele modelleringsprocessen. Det vil derfor være speci ceret, hvilken tidsrække der er brugt, når det er nødvendigt. En outlier er her en observation, hvis afstand fra et n dages moving average er større end
c standardafvigelser. I opgaven ses henholdsvis på n = 7 og n = 49 valgt som i Stoll and Wiebauer [2010], i øvrigt vælges c = 3. Observationer, der bryder med et af disse to ltre karakteriseres som outliers. Disse observationer erstattes med en lineær interpolation af de nærmeste punkter, der ikke selv er outliers. Filtret, som benyttes, er vist i ligning (4.1), hvor Ptm er et n dages moving average til tiden t.
15
(4.1)
Ft± = Ptm ± c · σ(Pt − Ptm )
Med outliers Uden outliers
25 20 10
15
Euro/MWh
30
35
40
I gur 4.1 ses gasprisen plottet både med og uden outliers.
2008
2009
2010
2011
2012
2013
2014
Årstal
Figur 4.1: Gasprisen illustreret både med og uden outliers.
4.2 Temperaturdata Temperaturdata er hentet fra Deutscher Wetterdienst og er frit tilgængeligt på http://www.dwd.de. Både gennemsnitstemperaturen samt den maksimale og minimale temperatur på dagen er indeholdt i data. I opgaven benyttes ligesom Stoll and Wiebauer [2010] temperaturen fra vejrstationen i Hannover. Dagstemperaturerne er tilgængelige helt tilbage til 1936, men vi har valgt at starte data den 1. oktober 1950 for at sikre datakvaliteten. Müller et al. [2015] og Stoll and Wiebauer [2010] benytter gennemsnitstemperaturen over dagen som input, mens Kremser and Rammerstorfer [2012] laver et gennemsnit af den maksimale og den minimale temperatur over dagen. Denne er givet ved ligning (4.2).
Tt =
M axT empt + M inT empt 2
16
(4.2)
I gur 4.2 er gennemsnittet af den maksimale og den minimale temperatur fra den 1. oktober
−20
−10
0
°C
10
20
30
1950 til den 31. marts 2014 illustreret.
1950
1960
1970
1980
1990
2000
2010
Årstal
Figur 4.2: Gennemsnittet af den maksimale og minimale dagtemperatur.
4.3 Oliedata I opgaven benyttes futuresprisen med kortest tid til udløb på en tønde ICE Brent råolie. Opbygningen af en futurekontrakt på olie adskiller sig fra en aktiefutures angående levering. Den sidste dag, en future kan handles, er 15 kalenderdage før den første kalenderdag i leveringsmåneden, hvis denne er en bankdag. Hvis det ikke er en bankdag benyttes foregående bankdag. Fremover vil olieprisen således betegne futuren med kortest tid til udløb. Olieprisen er hentet fra Bloomberg med koden CO1. Kontrakten handles kun på engelske bankdage, der er således ikke priser i weekenden eller på engelske helligdage. Olieprisen bliver brugt som input i modellerne via et gennemsnit af prisen over de sidste 180 dage. Dette gennemsnit består ca. af 130 bankdage. For at få daglige observationer lineært interpoleres prisen over weekender og helligdage. I gur 4.3 er olieprisen illustreret uden interpolation.
17
90 80 70 60 30
40
50
Euro
2008
2009
2010
2011
2012
2013
2014
Årstal
Figur 4.3: Futuresprisen med kortest tid til udløb på en tønde Brent råolie.
4.4 Gaslagerdata Gaslagerdata hentes fra Gas Infrastructure Europe (GIE). Opgivelsen af gaslagerdata er et relativ nyt tiltag. Fra september 2007 til september 2010 blev beholdningen opgivet på ugentlig basis. Fra januar 2010 startede opgivelsen af beholdningen på daglig basis. Vi ønsker at danne en tidsrække fra 2007 til 2014 med daglige observationer. For at have en konsistent tidsrække dannes en ugentlig tidsrække. Dvs. at der kun benyttes observationer fra mandage i det daglige data. Derefter er der foretaget en simpel lineær interpolation mellem de ugentlige observationer, hvilket resulterer i, at der er observationer for alle dage. I opgaven benyttes lagerdata opgivet i procent af den maksimale kapacitet. Det at benytte den procentvise lagerbeholdning gør, at der ikke tages højde for det nominelle niveau i gaslagrene. Dette har både positive og negative konsekvenser. På den negative side kan det nævnes, at en forøgelse i den maksimale kapacitet vil mindske risikoen for, at lagrene når et kritisk lavt niveau. Dette vil ikke blive fanget af procentsatsen. Grunden til, at niveauet ikke inddrages, er, at det ikke er klart fra data, (specielt omkring skiftet fra ugentligt til dagligt data), hvorvidt en positiv ændring i kapaciteten er et udtryk for en reel forøgelse i kapacitet, eller blot et udtryk for en forbedring i datakvaliteten. Man kunne forestille sig et lager, der først har opgivet data fra en given dato, men som hele tiden har været operationelt. Der er 18
således ikke tale om en reel forøgelse af lagerkapaciteten. Gaslagerdataet er plottet i gur
60 20
40
Procent
80
100
4.4.
2009
2010
2011
2012
2013
2014
Årstal
Figur 4.4: Gaslagerbeholdningen målt i procent af den maksimale kapacitet. Det ses, at lagerniveauet på normale år bunder omkring 40 %, mens at den i vinteren 20082009 og 2012-2013 er helt nede i 20%, før den vender. Faldet i 2009 skyldes en krise mellem Rusland og Ukraine, hvor faldet i 2013 var forårsaget af en lang og kold vinter.
Vi vil nu gå ind i opstillingen af vores spotprismodel.
19
5
Gasprismodel
Vi ønsker at opstille en model for gasprisen. I Stoll and Wiebauer [2010] opstilles en model, hvor temperaturen eksogent forklarer gasprisen. Denne model udvides i Müller et al. [2015] således, at olie inddrages som forklarende variabel. Vi vil her bygge oven på dette arbejde ved at inkludere gaslager beholdningen. I afsnit 5.1 vil modellen fra Stoll and Wiebauer [2010] kort gennemgås. Dette efterfølges af et afsnit, der omhandler Müllers udvidelse af modellen. Fra afsnit 5.3 og frem beskrives inddragelsen af gaslager informationen samt andre modi kationer, og selve modellen opstilles.
5.1 Stoll og Wiebauers gasprismodel Der vil være en øget efterspørgsel efter gas om vinteren, der skyldes opvarmning, hvilket giver højere gaspriser. Efterspørgslen efter gas er også højere på arbejdsdage end på weekendog helligdage. Dette skyldes, at virksomheder benytter mere energi på arbejdsdage. Sæsonudsvingene modelleres delvist deterministisk via en arbejdsdagsindikator og cosinus sinus funktioner, og delvist stokastisk via temperaturen.
Deterministisk sæson Den deterministiske sæson består af en indikator variabel s(t) = 1[t∈{Arbejdsdag}] , der antager værdien 1, hvis det er en arbejdsdag og ellers er 0. Udover denne indgår hhv. en cosinus og sinus funktion, der begge har en periode på et år. Tiden t måles i dage.
S(t) = β1 · 1[t∈{Arbejdsdag}] + β2 · cos
2π · t 365.25
+ β3 · sin
2π · t 365.25
(5.1)
Sæson via temperatur Både Müller et al. [2015] og Stoll and Wiebauer [2010] benytter normaliserede akkumulerede heating degree days til at beskrive temperaturens indvirkning på gasprisen. En heating degree day HHDt = max(15◦ − Tt , 0) angiver således, hvor mange grader temperaturen er
20
under 15◦ på dag t. Idéen er, at opvarmningen af bygninger påbegyndes, når temperaturen falder under 15◦ . Disse akkumuleres over vinteren. I en vinter med mange kolde dage vil der være en højere efterspørgsel efter gas i forhold til en vinter med færre koldere dage. Dette bør presse prisen på gas op. Vinteren starter den 1. oktober og løber de efterfølgende 181 dage. Lad HHDd,w angive heating degrees på dag d ∈ {1, ..., 182} i vinter nummer w. De akkumulerede heating degree days angives CHDDd,w , og er givet ved ligning (5.2).
CHDDd,w =
d X
HHDd,w
(5.2)
k=1
De akkumulerede heating degree days normaliseres ved at fratrække gennemsnittet af CHHDd,w fra de forrige år. Idéen er her, at en vinter, som er koldere end gennemsnittet, alt andet lige vil lede til højere gaspriser. Lad Λd,w angive de normaliserede CHHDd,w . w−1
Λd,w
1 X CHDDd,l = CHDDd,w − w − 1 l=1
Stoll and Wiebauer [2010] argumenterer for, at Λd,w kan bruges som en proxy for lagerbeholdningen i gaslagrene. En koldere vinter bør medføre lavere lagerbeholdninger som følge af højere efterspørgsel. Både Müller et al. [2015] og Stoll and Wiebauer [2010] lader i sommermånederne Λd,w falde lineært fra slutniveauet den 1. april tilbage til et niveau på 0 den 1. oktober. Müller et al. [2015] argumenterer for dette med en antagelse om, at gaslagrene fyldes hen over sommermånederne. Λd,w er illustreret i gur 5.1.
5.2 Müllers udvidelse med olie I Müller et al. [2015] gives to argumenter for at inkludere olie som forklarende variabel. En ren temperaturmodel fanger ikke, at efterspørgslen efter gas afhænger af de økonomiske konjunkturer. Dette blev specielt udstillet under krisen 2008-2009. Det første argument er således, at olie kan benyttes til at beskrive økonomiens tilstand. Det andet argument er, at lande som fx Tyskland importerer gas på kontrakter, der er indekseret på olieprisen. Værdiansættelsen af disse kontrakter sker ud fra olieformler. En olieformel indeholder tre 21
200 100 0 −100 −200
Normaliserede CHDD
−300 −400 2010
2011
2012
2013
2014
Årstal
Figur 5.1: Normaliserede akkumulerede heating degree days, der lineært interpoleres til 0 fra 1. april til 1. oktober. parametre. Averaging months angiver antallet af måneder, hvor der tages gennemsnit af olieprisen. Time lag er forskydningen mellem de måneder, der tages gennemsnit over, til de måneder, hvor olieformel-prisen er gældende. Validity months er antallet af måneder, hvor olieformel-prisen er gældende. I gur 5.2 ses et eksempel på olieformlen 3-1-1. Denne viser, at gennemsnittet over april, maj og juni er gældende for august.
| − April − | − Maj − | − Juni − | −Juli − | −August − | {z } | {z } | {z } | Averaging months
Time lag
Validity months
Figur 5.2: Olieformlen 3-1-1 Det er ikke umiddelbart oplagt, hvilken formel, der bedst beskriver gasprisen. Det ønskes at benytte olieformlen, der har den højeste forklaringsgrad. I markedet er de mest brugte formler 3-1-1, 3-1-3, 6-1-1, 6-1-3 og 6-3-3. Det er derfor oplagt at starte med at teste disse formlers forklaringsgrad, men da forskellige kontrakter prisfastsættes ud fra forskellige formler, kan den bedste lige så vel være en blanding af formlerne. Müller et al. [2015] viser, at det er formlen 6-0-1, der bedst beskriver gasprisen. Olieformlerne har den svaghed, at olieformel-prisen hopper. Det ses af gur 5.3, hvor olieformelprisen 6-0-1 hopper hver gang, der er gået en måned (antallet af validity months). Man kan med fordel udglatte olieformlerne. Dette vil hjælpe på olieformlernes forklaringsgrad af ga22
sprisen, da den ikke hopper. Desuden giver det også intuitiv mening i forhold til, hvordan markederne agerer. Olieformel-prisen er et gennemsnit af kendte priser. Derfor har markederne information om, hvad den næste olieformel-pris bliver, inden den træder i kraft. Markederne udglatter derfor selv overgangen fra en olieformel-pris til en anden. I stedet for at udglatte olieformlen 6-0-1, ses det på gur 5.3, at et 180 dages moving average er et godt
90
alternativ. Dette benyttes fremover.
60 40
50
Euro/bbl
70
80
6−0−1 180 days moving average
2008
2009
2010
2011
2012
2013
2014
Årstal
Figur 5.3: Olieformlen 6-0-1 og 180 dages moving average.
Opstilling af Müllers model Müllers model for gasprisen kan opstilles på formen i ligningen (5.3).
Pt = β0 +β1 ·1[t∈{Arbejdsdag}] +β2 ·cos
2π · t 2π · t +β3 ·sin +β4 Λt +β5 Ψt +XtG (5.3) 365.25 365.25
Pt er gasprisen, S(t) er den deterministiske sæsonvariabel, Λt er normaliserede akkumulerede heating degree days og Ψt 180 dages moving average. XtG er fejlledet, og denne stokastiske process modellerer Müller som en ARMA-process. Fjernes ledet β5 Ψt , er man tilbage ved modellen af Stoll and Wiebauer [2010]. Fittet af de to modeller er plottet i gur 5.4.
23
35 25 20 5
10
15
Euro/MWh
30
Gaspris Müller Stoll el al.
2009
2010
2011
2012
2013
2014
Årstal
Figur 5.4: Fittet fra modellerne holdt op mod gasprisen. Forklaringsgraden for modellen af Stoll et al. er 0.5229, mens Müllers model kommer op på 0.8397. Dette må siges at være en stor forbedring af modellens forklaringsgrad. Det er især faldet i krisen 08/09, som modellen af Stoll et al. har svært ved at forklare. Müllers model fanger derimod faldet i 08/09, og olieprisen kan derfor bruges som en proxy for økonomiens tilstand. Estimaterne for modellerne ses i tabel 13.1 i bilag.
5.3 Inddragelsen af gaslagerdata I udvidelsen af modellen indrages gaslagerbeholdningen procentvis ift. den maksimale kapacitet. Som med temperatur fratrækkes lagerbeholdningen den gennemsnitlige lagerbehold-
\ på samme dag de tidligere år. Idéen er, at det er afvigelser fra en normalbeholdning Storage ning, som forklarer udsving i gasprisen. Dette gør, at tidsrækken forkortes med et år, da man i det første år ikke har en normalbeholdning fra tidligere år som sammenligningsgrundlag. y−1
\ t= Storage
1 X Storaged,y y − 1 k=1
^ t = Storaget − Storage \ t Storage
24
Disse afvigelser forventes at have forskellig betydning alt efter tiden på året. En lavere beholdning end normalen forventes at have en større betydning for gasprisen om vinteren end om sommeren. For at fange dette inkluderes vekselvirkningen mellem lagerbeholdningen og cosinus og sinus fra den deterministiske sæson. Vekselvirkningerne tilfører en sæsoneekt af afvigelserne i gaslageret fra normal niveauet.
5.4 Temperaturchok Kremser and Rammerstorfer [2012] benytter temperaturdata på en anden måde end Müller et al. [2015]. Kremser and Rammerstorfer [2012] introducerer variablen temperaturchok. Et temperaturchok de neres som forskellen i degree days på den givne dag og gennemsnittet for den samme dag de 30 foregående år. Denne variabel beregnes hele året i modsætning til normaliserede akkumulerede heating degree days Λd,w (NHDD), der lineært interpoleres fra april til oktober. Derfor erstattes fodtegnet for vinter w med y for år. Desuden introduceres en cooling degree day. Ved høje gennemsnitstemperaturer nedkøles huse blandt andet ved brug af aircondition. Gas indgår som input i strømproduktionen til aircondition. For hver grad gennemsnitstemperaturen er over 21◦ stiger cooling degree målet med 1. Temperaturchok variablen er angivet i ligning (5.4). Bemærk, at heating degree days i Kremser and Rammerstorfer [2012] allerede starter ved 18◦ fremfor de 15◦ i Müller et al. [2015]. 30
1 X DDd,y−k Wt = DDd,y − 30 k=1
(5.4)
Degree days er summen af cooling degree days og heating degree days.
DDd,y = CDDd,y + HDDd,y = max{18◦ − T, 0} + max{T − 21◦ , 0}
(5.5)
I afsnit 5.5 benyttes temperaturchok fremfor NHDD. Dette skyldes, at den akkumulerede eekt ikke i samme grad er ønskværdig, da vi også inkluderer lagerdata. Desuden øges forklaringsgraden i regressionen, når temperaturchok benyttes. Grunden til dette kan være, at Müller et al. [2015] benytter NHDD som en proxy for lagerdata, hvilket er et direkte 25
input i modellen nedenfor.
5.5 Opstilling af model Nu opstilles en log-lineær model for gasprisen med det funktionelle udtryk givet i ligning (5.6).
2π · t 2π · t log Pt = β0 + β1 · 1[t∈{Arbejdsdag}] + β2 · cos + β3 · sin + β4 · log Ψt 365.25 365.25 2π · t 2π · t ^ + β8 · sin + XtG + β5 · Wt + Storaget β6 + β7 · cos 365.25 365.25 (5.6) Modellen skal benyttes til at simulere prisstier for gassen. For at sikre positive gasprissimulationer modelleres log til gasprisen. Tolkningen af koe cienterne bliver så, hvor mange procentpoint gasprisen ændres, når den forklarende variabel ændrer sig med 1 enhed. Dette indses ved at se på forholdet mellem gasprisen til to tidspunkter, hvor den eneste forskel er, at kovariat nummer i er 1 større på tid t2 end tid t1 .
Pt2 = eβi ≈ (1 + βi ) for βi tæt på 0 Pt1 Påvirkningen fra lageret er gjort sæsonafhængig ved at inkludere vekselvirkninger med cos og sin. I det venstre panel i gur 5.5 ses, hvor mange procentpoint gasprisen stiger, som følge af et fald i lagerchokket med et procentpoint. Grafen har således tiden t på 1. aksen og på 2. aksen haves:
2π · t 2π · t + βˆ8 · sin −1 exp −0.01 · βˆ6 + βˆ7 · cos 365.25 365.25
26
Lagerchok
10 0
Procentpoint
1.4 0.8
−20
1.0
−10
1.2
Procentpoint
1.6
20
1.8
Effekt på gasprisen
mar 13
jun 13
sep 13
dec 13
mar 14
sep 08
Dato
jun 09
mar 10
dec 10
sep 11
jun 12
mar 13
dec 13
Dato
Figur 5.5: Venstre panel angiver den tidsafhængige ændring af gasprisen i procentpoint, når lagerchokket falder med et procentpoint. Højre panel angiver det realiserede lagerchok og dennes eekt på gasprisen. Begge i procentpoint. Det ses, at hvis lagerbeholdningen falder med et procentpoint vil det i sommermånederne ikke have en ligeså stor eekt som i vintermåneder. Den laveste eekt ligger i starten af maj, hvor et fald med 1 procentpoint estimeres at forøge gasprisen med lidt over 0.8 procentpoint, hvor et tilsvarende fald vil forøge gasprisen med lidt over 1.8 procentpoint ved slut november start december. Olien indrages nu som log til de 180 dages moving-average. Hvis alle andre forklarende variable holdes fast, så vil relationen mellem Pt og Ψt være en potensfunktion. Derved vil en given procentuel stigning fra tid t1 til tid t2 på x procent, dvs. Ψt2 = Ψt1 (1 + x) give en procentvis stigning i gasprisen på Pt2 = Pt1 (1 + x)β4 . β -estimatets fortegn er positivt, hvilket intuitivt betyder, at hvis olieprisen stiger, så stiger gasprisen ligeledes. Denne tolkning stemmer nt overens med pointerne fra afsnit 5.2, hvor det blandt andet blev beskrevet, at Tyskland importerer gas på kontrakter, der er indekseret på olieprisen.
27
Intercept s(t) sin cos t Wt ^ t Storage log Ψt ^ t : cos Storage
Estimat -3.3 0.02 0.11 -0.08 -0.00 -0.01
Std. Error 0.07 0.01 0.00 0.00 0.00 0.00
t value -47.35 3.21 24.66 -17.67 -19.92 -12.36
Pr(>|t|) 0.00 0.00 0.00 0.00 0.00 0.00
-1.36 1.54
0.04 0.02
-31.33 82.74
0.00 0.00
-0.45
0.05
-8.65
0.00
-0.30
0.05
-6.50
0.00
^ t : sin Storage
Tabel 5.1: Parameterestimater for den lineære model af log-gasprisen. Modellen har en justeret R2 på 0.8978. Det ses i tabel 5.1, at fortegnet på β -estimatet for temperaturchok er negativt. En stigning i temperaturen vil medføre et fald i efterspørgslen. Dette resulterer i et fald i prisen. s(t) er som nævnt et 0-1 variabel, der beskriver, om det er hverdag. Denne påvirker gasprisen positivt, hvilket skyldes, at især virksomheder benytter mere gas på arbejdsdage end i weekender og
Gaspris Müller Model med lager
25 20 15 10 5
Euro/MWh
30
35
på helligdage. I gur 5.6 ses ttede værdier af vores model.
2009
2010
2011
2012
2013
Årstal
Figur 5.6: Fittede værdier af modellen i ligning (5.6)
28
2014
5.6 Modellering af XtG XtG er residualerne fra model (5.6). Her er det vigtigt at nævne, at der med residualerne menes forskellen mellem ttet fra modellen i ligning (5.6) og den faktiske gaspris, hvor outliers er indeholdt. Ser man på graferne for autokorrelationen og den partielle autokorrelation af residualerne, ses det, at det ikke er hvidstøj. Der ttes en AR(5) model dels på baggrund af det partielle korrelationsplot i gur 5.7, hvor det første insigni kante lag er lag 6. Endvidere bliver dette understøttet af, at en AR(5) har den laveste AIC. AIC
AR(1) -6537
AR(2) -6592
AR(3) -6658
AR(4) -6663
AR(5)
-6664
AR(6) -6663
AR(7) -6662
Tabel 5.2: Tabel over AIC for AR modeller. AR(5)-processen er opstillet i ligningen (5.7).
XtG
=
5 X
! G φk Xt−k
+ εt
(5.7)
k=1
εt er en hvidstøj process, og i tabel 5.3 ses estimaterne fra modellen. φ1 φ2 φ3 φ4 φ5
Estimat 0.71 0.02 0.14 0.03 0.05
Std. Error 0.02 0.03 0.03 0.03 0.02
Tabel 5.3: Estimaterne og standardafvigelserne på estimaterne. Det ses, at estimaterne for φ2 og φ4 er mindre en 2 standardafvigelser. De medtages grundet AIC og eekten på autokorrellationsplottene.
29
0
5
10
15
20
25
XG t
0.0 0.2 0.4 0.6 0.8
Partial ACF
0.4 0.0
ACF
0.8
XG t
30
0
5
10
Lag
20
25
30
Lag εt
0.0
0.00 −0.04
0.4
Partial ACF
0.8
0.04
εt
ACF
15
0
5
10
15
20
25
30
0
Lag
5
10
15
20
25
30
Lag
Figur 5.7: Første og andet panel angiver henholdsvis autokorrelationen og den partielle autokorrelation for residualerne (XtG ) i model (5.6). Tredje og fjerde panel er tilsvarende plot for εt I stil med Müller et al. [2015] og Stoll and Wiebauer [2010] estimeres en generaliseret hyperbolsk fordeling til støjen fra den autoregressive process. Klassen af generaliserede hyperbolske fordelinger er kontinuerte fordelinger på den reelle akse. De indeholder blandt andet t-fordelingen, variansgamma-fordelingen og den normale invers gaussiske fordeling. Disse tre fordelinger er karakteriseret ved, at de har tunge haler. AIC benyttes igen som kriterie for valg af fordeling. Det ses i tabel 5.4, at den symmetriske normale invers gaussiske fordeling (NIG) har den laveste AIC.
30
Model NIG ghyp t NIG ghyp t hyp hyp VG VG gauss
Symmetrisk sand sand sand falsk falsk falsk sand falsk sand falsk sand
lambda -0.50 -1.17 -1.73 -0.50 -1.20 -1.73 1.00 1.00 1.22 1.22
alpha.bar 0.68 0.53 0.00 0.68 0.51 0.00 0.40 0.40 0.00 0.00 Inf
mu 0 0 0 0 0 0 0 0 0 0 0
sigma 0.03 0.03 0.04 0.03 0.03 0.04 0.03 0.03 0.03 0.03 0.03
gamma 0 0 0 0 0 0 0 0 0 0 0
AIC -8166.40 -8166.16 -8165.77 -8164.68 -8164.38 -8164.08 -8151.26 -8149.59 -8145.91 -8144.29 -7819.42
Tabel 5.4: AIC vurdering af de generaliserede hyperbolske fordelinger.
31
llh 4086.20 4087.08 4085.89 4086.34 4087.19 4086.04 4078.63 4078.79 4075.95 4076.14 3911.71
6
Simulationsstudie
For at kunne generere prisstier for gasprisen, skal der opstilles modeller for temperatur, olie og gaslageret. Dette sker i afsnit 6.1-6.3. Simulationer fra disse modeller vil indgå som input til simulationen af gasprisen, der ses i afsnit 6.4.
6.1 Temperaturmodel Den gennemsnitlige dagstemperatur fra den 1. oktober 1950 til den 31. marts 2014 er illustreret i gur 4.2. Ikke overraskende ses det, at den har klare sæsonudsving med høje gennemsnitlige temperaturer om sommeren og lave om vinteren. For at modellere den deterministiske sæson benytter vi cosinus og sinus bølger med årlige frekvenser. Tiden t målt i dage starter den 1. oktober 1950.
Tt = β0 + β1 t + β2 cos
2πt 365.25
+ β3 sin
2πt 365.25
+ XtT emp
(6.1)
Koe cienterne i ligning (6.1) estimeres først med mindste kvadraters metode. Estimaterne ses i tabel 6.1. Residualerne XtT emp fra denne model er autokorreleret, hvilket ses i de to øverste paneler i gur 6.1. Det ses i det partielle autokorrelation plot, at der er 6 signi kante lags, hvorefter der ikke er den store partielle autokorrelation. Estimeres XtT emp som en AR(4) ses det i de nederste paneler af gur 6.1, at lag 5 og 6 ikke længere er signi kante. AR(4)-modellen er vist i ligning (6.2).
XtT emp
=
4 X
T emp φi · Xt−i + t
(6.2)
i=1
En samlet estimation af både sæson, trend og autoregressive komponenter kan ndes i tabel 6.2. De nederste to paneler i gur 6.1 viser, at der efter denne procedure ikke er nævneværdig autokorrelation tilbage. Et QQ-plot, baseret på, at støjen er normalfordelt, er vist i gur 6.2. Det ses, at halerne stikker lidt af fra linjen, men ikke nok til at lave en ny fordelingsantagelse.
32
βˆ0 βˆ1 βˆ2 βˆ3
Estimat 8.28 0.00 -8.08 3.04
Std. Error 0.05 0.00 0.03 0.03
t value 168.12 21.09 -232.10 87.31
Pr(>|t|) 0.00 0.00 0.00 0.00
Tabel 6.1: Estimater for modellen i ligning (6.1), hvor der ikke er taget højde for autokorrelationen i residualerne.
φˆ1 φˆ2 φˆ3 φˆ4 βˆ0 βˆ1 βˆ2 βˆ3
Estimat 0.922 -0.192 0.051 0.023 8.278 0.000 -8.081 3.040
Std. Error 0.007 0.009 0.009 0.007 0.146 0.000 0.103 0.103
0.6 0.4 0.2 0.0
Partial ACF
0.8
0.0 0.2 0.4 0.6 0.8 1.0
ACF
Tabel 6.2: Estimater for modellen i ligning (6.1), hvor der er taget højde for autokorrelationen i residualerne.
0
10
20
30
40
0
10
Series ARres
30
40
0.000 −0.010
Partial ACF
0.010
Lag
0.0 0.2 0.4 0.6 0.8 1.0
Lag
ACF
20
Series ARres
0
10
20
30
40
0
Lag
10
20
30
40
Lag
Figur 6.1: De øverste paneler er autokorrelationsplots for XtT emp . De nederste er for residualerne til AR(4) modellen. 33
Normal Q−Q Plot
0 −10
−5
Sample Quantiles
5
●●●●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ●●
●
●
●
−4
−2
0
2
4
Theoretical Quantiles
Figur 6.2: QQ-plot
6.1.1
Simulation af Temperatur
Modellen i ligning (6.1) simuleres ved hjælp af
R-funktionen simulate.Arima()
fra pakken
forecast af Hyndman. I gur 6.3 ses 100 simulationer illustreret ved den grå farve. Der er ligeledes plottet et gennemsnit af simulationerne, der viser formen på den deterministiske
−10
0
10
20
Realiseret temperatur Simulation Gennemsnit af simulationer
−20
Grader celsius
30
40
sæson.
2008
2010
2012
2014
Årstal
Figur 6.3: Simulerede stier for temperaturen.
34
I gur 6.4 er simulationerne af temperaturchok variablen illustreret. Disse skal bruges som input til gasprissimulationen og er dannet på baggrund af ovenstående temperatursimulation. Her er ligeledes plottet 100 simulationer. Gennemsnittet af simulationerne ligger og
Realiseret temperaturchok Simulation Gennemsnit af simulationer
−15
−10
−5
0
5
10
15
svinger omkring 0.
2008
2010
2012
2014
Årstal
Figur 6.4: Simulerede stier for temperaturchok variablen er beregnet med ligning (5.4).
6.2 Olie Mange har tidligere modelleret futurepriserne på olie; Schwartz [1997] samt Schwartz and Smith [2000] for at nævne nogle få. Da det kun er futuren med kortest tid til udløb, som indgår via et 180 dages moving average i modellen for gas, opstilles her bare en ARMA-GARCH model. I det øverste panel i gur 6.5 er
CO1 illustreret. Ved at foretage et Augmented Di-
ckey Fuller test kan nulhypotesen om, at tidsrækken har en enhedsrod, ikke forkastes, da testes p-værdi er 0.1583. Nulhypotesen, at tidsrækken har en enhedsrod og derved ikke er stationær, accepteres. Der tages log til tidsrækken, hvorefter den dierences. Tidsrækken kan nu tolkes, som afkastet på CO1 prisen. Lad αt angive afkastet mellem dag t − 1 og t.
CO1t = (1 + αt )CO1t−1 ∇ log CO1t = log(1 + αt ) ≈ αt 35
Den sidste approksimation er bedre des tættere αt er på nul. Autokorrelationsplottet og det
60 20
CO1
partielle autokorrelationsplot for ∇ log CO1t er angivet i de øverste to paneler af gur 6.6.
11000
12000
13000
14000
0.00 0.10 −0.15
∇logCO1
Årstal
2000
2002
2004
2006
2008
2010
Årstal
Figur 6.5: Tidsrækken og den log-dierensede tidsrække af CO1t Der modelleres en MA(1) model til ∇ log CO1t , da der i autokorrelationsplottet ses et signi kant lag. Autokorrelationen for MA(1) modellens residualer kan ses i de nederste to paneler af gur 6.6. Modellen skal ligeledes kontrolleres for volatilitetsclustering i tidsrækken. Kontrollen udføres ved at kigge på autokorrelationen og den partielle autokorrelation i de kvadrerede residualer fra MA(1)-modellen. Disse ses i gur 6.7 og viser, at der er volatilitetsclustering i data. Dette modelleres med en Garch(1,1) på baggrund af BIC. Se tabellerne 6.3 og 6.4. Estimationen er foretaget med R-pakken fGarch af Wuertz et al. [2013]. Den samlede model kan skrives som:
∇ log CO1t = µ + θυt−1 + υt υt = ε t σ t
(6.3)
2 σt2 = ω + α1 ε2t−1 + β1 σt−1
I første omgang var fordelingsantagelsen for εt ∼ N (0, σt2 ). Dette gav anledning til QQplottet i panelet nederst til venstre i gur 6.7. Her ligger punkterne ikke helt pænt på linjen. 36
Specielt ikke i den nedre hale. Derved ændres fordelingsantagelsen til, at εt følger en skæv tfordeling. QQ-plottet med denne antagelse er vist nederst til højre i gur 6.7, hvor punkterne placerer sig bedre på linjen. Jævnfør bilag tabel 13.2 for at se estimaterne for den skæve t-fordeling.
i /j
1 2 3 4
0 -4.878 -4.914 -4.937 -4.965
1
-5.037
-5.035 -5.032 -5.030
2 -5.035 -5.033 -5.031 -5.029
3 -5.034 -5.034 -5.033 -5.030
4 -5.032 -5.031 -5.029 -5.028
Tabel 6.3: BIC for Garch(i,j ) speci kationer, hvor εt er antaget normalfordelt.
i /j
1 2 3 4
0 -4.953 -4.971 -4.983 -4.997
1
-5.057
-5.055 -5.052 -5.050
2 -5.055 -5.053 -5.050 -5.048
3 -5.053 -5.054 -5.052 -5.049
4 -5.051 -5.050 -5.048 -5.046
−0.04
0.00
Partial ACF
0.4 0.0
ACF
0.8
0.04
Tabel 6.4: BIC for Garch(i,j ) speci kationer, hvor εt er antaget at følge en skæv t-fordeling.
0
5
10
15
20
25
30
35
0
5
10
Series res
20
25
30
35
25
30
35
Lag
0.0
0.00 −0.04
0.4
Partial ACF
0.8
0.04
Lag
ACF
15
Series res
0
5
10
15
20
25
30
35
0
Lag
5
10
15
20 Lag
Figur 6.6: ACF og PACF for ∇ log(CO1 ) ses i de to øverste paneler. I de to nederst er ACF og PACF plottet for residualerne til MA(1)-modellen. 37
0.20 0.00 10
15
20
25
30
35
0
20
25
Lag
−2
0
2
4
● ●
−6
Theoretical Quantiles
35
●
2
Sample Quantiles 4
30
●●●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●●● ● ●
0
4 2 0
15
Skeewed t Q−Q Plot
●
−4
10
Lag ●
−6
5
Normal Q−Q Plot
−6 −4 −2
5
●●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●
−6 −4 −2
Sample Quantiles
0.10
Partial ACF
0.8 0.4 0.0
ACF
0
−4
−2
0
2
4
Theoretical Quantiles
Figur 6.7: I de to øverste paneler er ACF og PACF for de kvadrerede residualer vt2 plottet. I de nederste paneler er QQ-plottet for residualerne ved en normalfordeling og en skæv t-fordeling.
6.2.1
Simulation af oliepriser
Modellen i ligningssystemet (6.3) simuleres ved hjælp af R-funktionen garchSim() fra pakken
fGarch af Wuertz et al. [2013]. Der foretaget 100 simulationer. Disse be nder sig i spændet mellem 15 og 150. I gur 6.8 ses, at gennemsnittet for simulationerne har en nedadgående tendens. I gur 6.9 ses 100 simulationer af det 180 dages moving averages. Det bemærkes, at det ikke er de samme 100 simulationer som ovenfor. Dette ses tydeligst via simulationen, der bevæger sig over 160 EURO. Nedenstående simulation benyttes som input til gasprissimulationen.
38
160 140 100 80 20
40
60
Euro
120
Realiseret oliepris Simulation Gennemsnit af simulationer
2008
2010
2012
2014
Årstal
Realiseret 180 dages moving average Simulation Gennemsnit af simulationer
100 80 20
40
60
Euro
120
140
160
Figur 6.8: Simulerede stier for olieprisen.
2008
2010
2012
2014
Årstal
Figur 6.9: Simulerede stier for 180 dages moving average af olieprisen.
6.3 Gaslager model Gasbeholdningen, er som tidligere nævnt, målt som den procentvise beholdning af den maksimale kapacitet. Den ideelle model for lagerbeholdningen ville blandt andet afhænge af gasudbuddet, primært fra Rusland og Norge. Ved et lavt udbud, vil lagrene gå ind og dække efterspørgslen og det vil medføre et fald i lagerbeholdningen. Det ses på gur 4.4, at 39
lagerbeholdningen var meget lav i starten af 2009 og 2013. Det lave niveau i 2009 skyldtes, at Rusland den 7. januar 2009 stoppede sin eksport af gas til Europa. Det lave niveau i 2013 tilskrives en lav import af LNG, hvilket gjorde, at en større del af energiforbruget skulle dækkes med gas. Det er om noget endnu sværere at kvanti cere og modellere disse faktorer end gasbeholdningen selv. Derfor vil gasbeholdningen her udelukkende modelleres ved hjælp af deterministiske sæsonkomponenter. De deterministiske sæsonkomponenter er cosinus og sinus bølger både med helårlig og halvårlig frekvens. Argumentet for at bølger med halvårlige frekvenser benyttes er, at det fysisk tager længere tid at fylde et lager end at tømme det. Den halvårlige frekvens tilfører denne asymmetri til svingningerne. I modelleringen af gaslagerbeholdningen transformeres procentdelen af den maksimale kapacitet med logit. Transformationen kan således antage alle værdier på den reelle akse. Dette gøres, da man så kan lave en lineær model, hvor man har sikret sig, at gaslagerbeholdningen forbliver i intervallet [0,1], når der transformeres tilbage. Modellen estimeres ud fra ugentlige observationer ved hjælp af
R-funktionen Arima()
fra pakken
forecast af Hyndman. Først estimeres det lineære t,
hvorefter ARMA(2,1) modellen estimeres på residualene fra 29. marts 2010 og frem. ARMA(2,1) er valgt grundet BIC i tabel 6.5 og ACF og PACF vist i gur 6.10. Modellen burde estimeres på alle residualerne tilbage fra 2007. Dette er ikke valgt, da volatiliteten falder med årene. Det kan skyldes, at det nominelle lagerniveau løbende stiger samt forbedringen af datakvaliteten, der blev beskrevet i afsnit 4.4. Modellen er angivet i ligning (6.4), hvor
t ∼ N (0, σ 2 ) og parameterestimaterne er givet i tabel 6.6. log
Storaget 1 − Storaget
2 X 2kπ · t 2kπ · t + β2k sin + XtLager = β0 + β2k−1 cos 365.25 365.25 k=1
Lager Lager XtLager = φ1 Xt−1 + φ2 Xt−2 + θ1 εt−1 + εt
(6.4)
40
AR(i )\MA(j ) 1 2 3 4 5
1 -633.38
-659.35
-654.31 -649.48 -645.39
2 -650.65 -654.70 -649.39 -647.09 -641.84
3 -655.86 -650.60 -646.40 -640.78 -642.50
4 -650.56 -645.98 -641.93 -636.63 -632.99
5 -646.12 -641.54 -639.30 -633.08 -628.10
Tabel 6.5: BIC for ARMA(i,j ) modeller.
φˆ1 φˆ2 θˆ1 βˆ0 βˆ1 βˆ2 βˆ3 βˆ4
Estimat 1.58 -0.60 -0.22 0.90 0.62 1.40 0.20 0.09
Std. Error 0.08 0.08 0.09 0.03 0.04 0.04 0.04 0.10
0.5 −0.5
0.0
Partial ACF
0.4 0.0
ACF
0.8
1.0
Tabel 6.6: Parameterestimater til modellen for gaslager givet ved ligning (6.4).
0
5
10
15
20
25
5
Series resid(ARMA21)
15
20
25
Lag
0.00 −0.10
0.0
0.4
Partial ACF
0.8
0.10
Lag
ACF
10
Series resid(ARMA21)
0
5
10
15
20
25
5
Lag
10
15
20
25
Lag
Figur 6.10: Autokorrelationplots og partielle autokorrelationsplots. De øverste paneler er af Xt , dvs. før en ARMA model er ttet. I de nederste to paneler er det af εt , dvs. efter en ARMA(2,1) er ttet. 41
6.3.1
Simulation af gaslagerbeholdning
Modellen i ligning (6.4) simuleres ved hjælp af pakken
R-funktionen simulate.Arima()
også fra
forecast af Hyndman. Dette giver ugentlige observationer på logitskalaen. Disse
observationer transformeres tilbage, således at de giver en værdi i procent. For at få daglige observationer lineært interpoleres mellem de ugentlige procentværdier. Simulationerne er vist
60 20
40
Procent
80
100
i gur 6.11. Det ses, at simulationernes amplitude ikke er den samme for alle simulationerne.
2008
2010
2012
2014
Årstal
Figur 6.11: Simulerede stier for gaslagerbeholdningen.
6.4 Simulation af gaspriser Alle de nødvendige input til modellen i ligning (5.6) er nu simuleret. Til at simulere gasprisen benyttes
R-funktionen Predict()
fra pakken
stats af R Core Team [2014]. Desuden
simuleres fejlledet, der skal tillægges prædiktionen, ved at benyttet funktion
fejl_sim().
Denne funktion er selvlavet og kan ndes i bilag under afsnit 13.4. Sidst benyttes eksponentialfunktionen til at transformere data fra log(Pt ) til Pt . Simulationen af gasprisen ses i gur 6.12.
42
50 30 20 0
10
Euro/MWh
40
Realiseret pris Simulation Gennemsnit af simulationer
2008
2010
2012
2014
Årstal
Figur 6.12: Simulerede stier af gaspriser. Gennemsnittet af gasprissimulationerne indikerer, at gasprisen vil ,starte med at falde. Dette sker primært som en konsekvens af, at de oliepriser, der er tættest på at falde ud af vinduet i det 180 dages moving average, er højere end olieprisen ved indgangen til april, som alt andet lige må give en indikation af olieprisen i den nære fremtid. Således indikerer faldet i olieprisen, at gasprisen vil falde i sommermånederne, hvorefter gennemsnittet igen stiger i vintermånederne. Dette afrunder afsnittet om modellerne, der simuleres på baggrund af. Nogle af modellens kritiske antagelser og egenskaber vil nu blive diskuteret, inden simulationerne kalibreres til forward kurven og fokus rettes mod prisfastsættelsen af gaslager kontrakter.
43
7
Backtesting og diskussion af model
I dette afsnit vil vi undersøge, hvorledes modellen præsterer out of sample. Dette gør vi i afsnittet 7.1. Dernæst vil vi i afsnit 7.2 belyse problemstillingen vedrørende kausalitet mellem gaslagerniveauet og gasprisen. Afslutningsvis vil vi i afsnit 7.3 kritisk diskutere de simulerede gaspriser, der leder til et forslag om en reduktion af de simulerede stier.
7.1 Backtesting I dette afsnit testes spotmodellen out of sample. Dette sker ved, at modellen estimeres ud fra en delperiode. Estimationen bliver så brugt til at forklare den delmængde af data, som ikke er inkluderet i estimationen. Det komplette datasæt løber fra perioden den 17. september 2008 til den 31. marts 2014. Datasættet inddeler vi i følgende perioder: Periode Periode Periode Periode Periode
1 2 3 4 5
Delperiode oversigt 17. September 2008 til 31. 1. April 2010 til 31. 1. April 2011 til 31. 1. April 2012 til 31. 1. April 2013 til 31.
Marts Marts Marts Marts Marts
2010 2011 2012 2013 2014
Tabel 7.1 Et mål for, hvor godt modellen præsterer out of sample, er Root Mean Squared Error RMSE, der er givet i ligning (7.1). Dette er kvadratroden af den gennemsnitlige kvadrede fejl MSE. MSE er givet ved summen af variansen på estimatet samt estimatets bias sat i anden: MSE(Pˆt ) = Var(θˆt ) + Bias(θˆt )2 . For et middelret estimat vil RMSE altså være sammenligneligt med standardafvigelsen for residualerne σ ˆ in sample.
RMSE =
√
v u n uX (θˆt − θt )2 MSE = t n t=1
(7.1)
Den udvidede gaspris model givet i ligning (5.6) estimeres nu på en række forskellige delperioder af tidsrækken og efterprøves out of sample. Den estimerede gas pris er beregnet ved
44
ligning (7.2).
2π · t 2π · t ˆ ˆ ˆ ˆ + β3 · sin + βˆ4 · log Ψt log Pˆt = β0 + β1 · 1[t∈{Arbejdsdag}] + β2 · cos 365.25 365.25 2π · t 2π · t ^ t βˆ6 + βˆ7 · cos + βˆ5 · Wt + Storage + βˆ8 · sin 365.25 365.25 (7.2) I første out of sample test estimeres modellen på periode 1 og afprøves på periode 2, dernæst estimeres modellen på periode 1 og 2 og testes på periode 3 og så fremdeles. I tabel 7.2 er estimatet for residualernes standardafvigelsen og RMSE i out of sample perioden givet. Det ses, at modellen præsterer dårligt med en RMSE på 13.7596, når den kun estimeres på baggrund af periode 1. Den dårlige præstation kan skyldes en række af ting. Eksempelvis indgår en række årlige sæsonkomponenter i modellen i form af cosinus- og sinussvingningerne, som her estimeres på halvandet års data. Man kan formode, at estimaterne tilpasser sig for meget til estimationsperioden, dvs. modellen over tter. Det ses i tabellen, at RMSE falder når estimationsperioden øges. En gra sk illustration af prædiktionerne i out of sample perioderne mod de realiserede værdier ses i gur 7.1. Estimationsperiode 1 1-2 1-3 1-4
Out of sample 2 3 4 5
σ ˆ (log(Pˆt )) RMSE(log(Pˆt )) σ ˆ (Pˆt ) RMSE(Pˆt ) 0.1020 1.1642 1.6277 13.7596 0.1042 0.3035 1.8277 8.4498 0.1193 0.1409 2.1493 3.9356 0.1143 0.1181 2.0600 2.8813
Tabel 7.2: Out of sample præstation for modellen, når estimationsperioden løbende forlænges. Grundet den mindre gode præstation fra modellen ved en lille estimation periode ses der nu på, hvor godt modellen præsterer, når den estimeres på 4 af de 5 delperioder og afprøves på den sidste. RMSE er angivet i tabel 7.3 og out of sample perioderne er illustreret i 7.2, hvor man kan se den forbedrede RMSE.
45
Estimationsperiode 2-5 1,3-5 1-2,4-5 1-3,5
Out of sample 1 2 3 4
σ ˆ (log(Pˆt )) RMSE(log(Pˆt )) σ ˆ (Pˆt ) RMSE(Pˆt ) 0.0635 0.2928 1.3208 3.9104 0.1016 0.1528 1.8739 2.8638 0.1144 0.0889 2.0217 2.0617 0.1136 0.1046 1.9982 2.8463
Tabel 7.3: Out of sample præstation, hvor estimationsperioden består af 4 af delperioderne og out of sample perioden består af en enkelt delperiode. Det ses, at når der benyttes længere delperioder i estimationen præsterer modellen generelt bedre. RMSE er størst i delperioden 2008-2010. Dette er en periode med store udsving, hvor gasprisen falder fra et niveau på over 30 EUR/MWh til under 10 EUR/MWh. Tilsvarende er
Jul 12
Okt 12
35 25 15 35
Jul 11
Okt 11
Jan 12
Okt 13
Jan 14
Gaspris i Eur/MWh
30
35 30 25 Apr 12
Apr 11
25
Jan 11
45
dat01.04.2011_01.04.2012[, 1] Okt 10
20
Jul 10
Gaspris i Eur/MWh
Gaspris i Eur/MWh
15
Apr 10
dat01.04.2013_01.04.2014[, 1]
30 5
10
15
20
25
Gaspris i Eur/MWh
20
dat01.04.2012_01.04.2013[, 1]
dat01.04.2010_01.04.2011[, 1]
RMSE mindst i delperioden 2011-2012, hvor gasprisen svinger mellem 20 og 26 EUR/MWh.
Jan 13
Apr 13
Jul 13
Figur 7.1: Out of sample præstation for modellen, hvor estimationsperioden baserer sig på alle delperioder før out of sample perioden.
46
26 24
Gaspris i Eur/MWh
20 18 16
5
14
10
15
20
25
Gaspris i Eur/MWh
22
30
35
Gaspris i Eur/MWh
Dec 08
Mar 09
Jun 09
Sep 09
Dec 09
Mar 10
Apr 10
Gaspris i Eur/MWh
Jul 10
Okt 10
Jan 11
Okt 12
Jan 13
Gaspris i Eur/MWh
32 30 28 22
18
24
26
Gaspris i Eur/MWh
26 24 22 20
Gaspris i Eur/MWh
28
34
30
Sep 08
Apr 11
Jul 11
Okt 11
Jan 12
Apr 12
Jul 12
Figur 7.2: Out of sample præstation, hvor alle andre delperioder end den i grafen er benyttet i estimationen.
7.1.1
Out of sample på lageråret 14/15
Modellen fra ligning (5.6) estimeres på perioden fra 17. September 2008 til den 31. Marts 2014 og ttes til lageråret 14/15. RMSE er her 4.1797 EUR mod 3.9006 EUR in sample. På log skalaen fås RMSE på 0.2237, hvilket skal holdes op mod de 0.1091 in sample. RMSE på EUR skalaen stiger således i forhold til de 2.8813, der blev opnået i lageråret 13/14. Man ville ellers forvente, at RMSE falder, når estimationsperioden forlænges. Hvorfor er det ikke tilfældet her? Modellens øgede fejl må primært tilskrives olieprisens styrtdyk, jævnfør guren i bilag 13.1. Styrtdykket er forårsaget af oliekartellet OPECs priskrig mod skiferolie. Udvindingen af skiferolien i USA er ikke rentabel, hvis olieprisen holdes på et lavt niveau. OPEC sikrer et lavt prisniveau ved skabe et øget udbud. Men tager modellen ikke højde for det? Som tidligere nævnt angives det i Müller et al. [2015], at olieprisen både påvirker gasprisen via indeksering, men at den også fungerer som en konjunkturvariabel.
47
26 20 14
16
18
Eur/MWh
22
24
Realiseret Gaspris Estimeret Gaspris Forwardkurven d. 31. marts 2014
Maj
Jul
Sep
Nov
Jan
Mar
2014−2015
Figur 7.3: Out of sample performance i lageråret 2014-2015, hvor den udglattede forwardkurve også er angivet. I starten af estimationsperioden tager olieprisen et drastisk fald i forbindelse med nanskrisen. Dog er konjunkturerne vidt forskellige ved faldet i olieprisen forårsaget af krisen og faldet forårsaget af priskrigen mod skiferolie. Da modellen ikke på andre måder bliver kontrolleret for konjunkturer, vil den reagere ens på et givent olieprisfald i begge tilfælde på trods af forskellen i konjunkturer. Dette gør, at modellen overvurderer olieprisfaldets eekt på gasprisen. Modellen forventer således, at gasprisen vil falde mere end, hvad den rent faktisk gjorde. Denne information er dog ikke kendt, når modellen opstilles den 31. marts 2014. Derfor er denne information ikke benyttet til at njustere modellen. Det danner mere et grundlag for at efterrationalisere modellen med henblik på at forbedre modellen til en prisfastsættelse af en gaslagerkontrakt i gasåret 15/16.
48
7.2 Kausalitet Et problem ved at benytte lagerbeholdningen som en eksogen forklarende variabel er, at der kan være feedback-mekanismer mellem denne og gasprisen. Hvis lagerbeholdningen er lav, vil der alt andet lige være et lavt udbud, hvilket vil resultere i højere priser. Omvendt hvis prisen er høj af en anden årsag, da vil man som ejer af lageret kunne pro tere på de høje priser ved at sælge gas. Gassalget vil så sænke beholdningen. På baggrund af dette forventes kausaliteten ikke umiddelbart at være entydig, hvilket gør at estimaterne i tabel 5.1 kan være biased.
En metode til at få unbiased estimater er ved at opstille en simultan ligningsmodel for gasprisen og lagerbeholdningen. Denne kan estimeres ved two stage least squares. De endogene variable er gasprisen og gaslagerbeholdningen. Olie, temperatur og arbejdsdagsindikatoren fungerer som instrumentvariable for lagerbeholdningen. I ligningen for lagerbeholdningen er det kun cosinus- og sinusbølgerne med en halvårlig cyklisk frekvens, som ikke indgår i ligningen for gasprisen. Disse benyttes som instrumentvariable for gasprisen. De estimerede parametre er angivet i tabel 7.4 og en graf af modellens t for gasprisen er vist i gur 7.4.
log Pt = β0 + β1 1[t∈{Arbejdsdag}] + β2 cos
2π · t 365.25
\ t) + + β5 Wt + β6 (Storaget − Storage log
Storaget 1 − Storaget
+ β3 sin
2π · t 365.25
+ β4 log Ψt
(7.3)
XtG
2 X 2kπ · t 2kπ · t = β0 + β2k−1 cos + β2k sin +β5 log Pt +XtLager 365.25 365.25 k=1 (7.4)
49
Ligning Gaspris Gaspris Gaspris Gaspris Gaspris Gaspris Gaspris Gaspris Lager Lager Lager Lager Lager Lager
Eekt Intercept sin1 cos1 s(t) t log Ψt Wt ^ Storage Intercept sin1 cos1 sin2 cos2 log Pt
Estimat -3.17 0.08 -0.09 0.02 -0.00 1.49 -0.01
Std. Error 0.09 0.01 0.01 0.01 0.00 0.02 0.00
t værdi -33.75 14.60 -16.99 3.18 -12.23 60.25 -11.60
Pr(>|t|) 0.00 0.00 0.00 0.00 0.00 0.00 0.00
-0.95 1.15 0.70 1.31 0.20 -0.05 -0.11
0.09 0.09 0.01 0.01 0.01 0.01 0.03
-10.88 12.47 51.60 97.85 15.19 -3.94 -3.70
0.00 0.00 0.00 0.00 0.00 0.00 0.00
Tabel 7.4: . Residualerne XtLager fra ligning (7.4) er illustreret i gur 7.5. Denne process ser ikke stationær ud. Det blev forsøgt at dierense tidsrækken. Dette gav en tidsrække, der oscillerede omkring nul. Det skabte dog det problem, at simulationerne af XtLager divergerede og ikke var stationære. Dette medførte, at lagerbeholdningens simulationer til tider forblev 100% eller 0%. Dette er ikke ønskværdigt og på baggrund af dette vil modellen, hvor gaslagerbeholdningen eksogent beskriver gasprisen, benyttes til simulationsstudiet. Desuden er der ikke markante forskelle på 2SLS- ttet i gur 7.4 og ttet af den biased model fra ligning (5.6) illustreret i gur 5.6. Dette kan indikere, at biaset ikke er så stort.
50
35 20 5
10
15
Euro/MWh
25
30
Gaspris 2SLS−fit
2009
2010
2011
2012
2013
2014
Årstal
1.0 0.5 0.0 −0.5 −1.0
Logit(Storage) residualer
1.5
Figur 7.4: Fittet af modellen estimeret med 2 stage least squares.
2009
2010
2011
2012 Årstal
Figur 7.5: Residualerne XtStorage
51
2013
2014
7.3 Reduktion i antallet af stier grundet ikke stationær oliemodel I afsnit 6 er der foretaget et simulationsstudie. Der er simuleret 30.000 stier for modellen, som kan benyttes til prissætningen af et lager. Prisfastsættelsen ved Least Squares Monte Carlo bliver bestemt ud fra de simulerede stier. Derved afhænger prisfastsættelsen af simulationsmodellen, hvilket medfører, at værdiansættelsen er behæftet med modelrisiko. Derfor er det vigtigt, at man vurderer, at modellen genererer troværdige simulationer. Desværre er der for vores model en lille mængde af stierne, der mangedobler gasprisen over året. Ligeledes er der stier, som går mod nul. Prisen den 1. april 2014 ligger omkring 20 EUR/MWh og fra 2008 og frem har prisen ikke været over 40 EUR. Gennemsnittet for de 2.5% højeste priser på sluttidspunktet er på 117.5352 EUR/MWh, hvilket er ca. seks gange højere end prisen den 1. april 2014. Grundene til de ekstreme gaspriser kan være mange. Modellen for gasprisen er opstillet i ligning (5.6) og her ses det, at gasprisen blandt andet afhænger af olieprisen, temperaturen og lagerniveauet. En ekstrem observation kan forekomme, hvis mindst en af disse simulationsmodeller giver ekstreme værdier. Det sker primært som en konsekvens af den ikke stationære oliemodel, der er modelleret med en enhedsrod. Der ligger således en ekstrem olie sti bag stort set alle de ekstreme gaspriser, som enten mangedobles eller går mod nul. Det er vores vurdering, at en forbedring af den samlede simulation kan opnås ved at opstille en stationær oliemodel, der fx konvergerer mod den næste månedlige kontrakt. Hvis gasprisen reelt skulle blive meget høj, ville man også kunne substituere gas med fx kul i strømproduktionen.
De ekstreme gaspriser kan medfører en misvisende værdiansættelse af lageret. Derfor har vi besluttet at lave en efterbehandling af simulationerne, hvor vi fjerner de simulationer, der bliver urealistisk høje og lave. Vi vil i afsnit 11 se prisfastsættelsen henholdsvis med og uden reduktion i stierne. Reduktionen af stierne er foretaget med udgangspunkt i den maksimale gaspris for hver sti. De 2.5% største og laveste observationer ndes, og de tilhørende stier bliver derefter fjernet. Det svarer til at fjerne 1.500 af de 30.000 simulerede stier.
52
8
Kalibrering af spotsimulationer til forwardpriserne
I dette afsnit vil vi kalibrere simulationerne til den observerede forwardkurve. Forwardkurven for gas består af kontrakter, der har forskellige leveringsperioder. Prisen på disse kontrakter kan ses som en vægtet sum af nutidsværdien af priserne over leveringsperioden. I afsnit 8.1 vil vi udglatte den observerede forwardkurve. Dette sker ud fra en betragtning om, at dagsprisen på gas ved skift mellem to perioder må ligge tæt på hinanden. Udglatningen foretages således for at undgå spring i simulationerne ved overgangen mellem to leveringsperioder. I afsnit 8.2 modi ceres simulationerne til den observerede forwardkurve. Dette sker på den sidste dag inden lagerårets start den 1. april, hvilket vil medføre, at gennemsnittet for simulationerne på en given dag vil være lig den udglattede forwardpris.
8.1 Udglatning af forwardkurven Lad F (t0 , T b , T e ) være den observerede forwardpris på tidspunkt t0 for forwardkontrakten, der leverer gas mellem tidspunkterne T b og T e . Forwardkontrakten leverer gassen stabilt over leveringsperioden [T b , T e ]. Køber man forwardkontrakten til tid t0 og følger handelsstrategien, hvor man sælger gassen igen lige så snart, at man har fået den leveret, da fås payoet givet i (8.1). Dette kan ses som payoet for en swapkontrakt, hvor forwardprisen
F (t0 , T b , T e ) udgør det faste ben. e
T X ti
(8.1)
Pti − F (t0 , T b , T e )
=T b
Hvis kontrakten afregnes ved afslutningen af leveringsperioden vil prisen kunne ndes ud fra relationen i ligning (8.2).
EQ
e
T X ti
Pti − (T e − T b )F (t0 , T b , T e )|Ft0 = 0
=T b
53
(8.2)
I denne isoleres den observerede forwardpris. e
T X 1 EQ [Pti |Ft0 ] F (t0 , T , T ) = e T − Tb b b
e
(8.3)
ti =T
Det ønskes at udglatte den observerede forwardkurve, der er ad over leveringsperioden. Se gur 8.1. Dette gøres ved at nde en glat og kontinuert funktion for forward prisen f (t0 , t), der fungerer som estimat for prisen på gas, der købes på tid t0 og leverer gassen på tid t.
f (t0 , t) skal opfylde relationen i ligning (8.4). b
Z
e
Te
F (t0 , T , T ) = Tb
1 f (t0 , t)dt Te − Tb
(8.4)
Således kalibreres de simulerede spotpriser til forwardkurven, hvis de opfylder relationen i ligning (8.5).
Z
t+1
f (t0 , s)ds = EQ [P (t)|Ft0 ]
(8.5)
t
Vi har benyttet approksimationen
R t+1 t
f (t0 , s)ds ≈ f (t0 , t). Vi har udglattet forwardkurven
ved brug af metoden beskrevet i Benth et al. [2007], hvor den funktionelle form af f (t0 , t) er givet ved summen af fjerdeordens spline. Denne er givet ved ligning (8.6) og et sæson element, hvor n angiver antallet af perioder med sit eget funktionelle udtryk. Analysen i Benth et al. [2007] viser, at den apriorisk valgte sæsonkomponent har en negligerbar betydning for udglatningen af forwardkurven. Derfor har vi ikke inkluderet en sæsonkomponent, men nøjes med fjerdeordens splinen. Udglatningen er illustreret i gur 8.1.
(t) =
a1 t4 + b1 t3 + c1 t2 + d1 t + e1 , a t4 + b t3 + c t2 + d t + e , 2
2
2
2
2
.. . an t4 + bn t3 + cn t2 + dn t + en ,
for t ∈ [T0 , T1 ] for t ∈ [T1 , T2 ]
for t ∈ [Tn−1 , Tn ]
x = [a1 b1 c1 d1 e1 a2 b2 c2 d2 e2 · · · an bn cn dn en ]0 54
(8.6)
(8.7)
Koe cienterne, der er indeholdt i x, ndes ud fra et maksimalt udglatningskriterie. Mere præcist løses minimeringsproblemet (8.8) begrænset af bibetingelserne (8.9)-(8.14):
Z
Tn 00
[ (t; x)] dt = min
min
x
T0
n Z X j=1
Tj
[12aj t2 + 6bj t + 2cj ]2 dt
(8.8)
Tj−1
23 22
Eur/MWh
24
25
x
2
●
21
Observed Smoothed Maj
Jul
Sep
Nov
Jan
Mar
Figur 8.1: En udglatning af forwardkurven for lageråret 2014-2015. Udglatningen er baseret på de observerede forwardpriser den 31. marts 2014.
(aj+1 − aj )t4j + (bj+1 − bj )t3j + (cj+1 − cj )t2j + (dj+1 − dj )tj + (ej+1 − ej ) = 0
(8.9)
4(aj+1 − aj )t3j + 3(bj+1 − bj )t2j + 2(cj+1 − cj )tj + (dj+1 − dj ) = 0
(8.10)
12(aj+1 − aj )t2j + 6(bj+1 − bj )tj + 2(cj+1 − cj ) = 0
(8.11)
0 (tn , x) = 0 Z Tie (t)dt = FiC
(8.12) (8.13)
e1 = P 0
(8.14)
Tib
Bibetingelserne (8.9), (8.10) og (8.11) sikrer, at værdien af (t), samt dens første og anden 55
a edte er ens fra højre og venstre i overgangene mellem de forskellige kontraktperioder. Derved giver (8.9)-(8.11) anledning til 3(n − 1) bibetingelser. Bibetingelsen (8.12) sikrer, at der ved udgangen af den sidste periode er en vandret tangenthældning. Bibetingelse (8.13) sørger for, at den kontinuerte forwardfunktion integreret over en given leveringsperiode er lig med de observerede forwardpriser. Denne bidrager med m bibetingelser, hvor m er antallet af forwardkontrakter. Den sidste bibetingelse (8.14) har vi tilføjet til rækken af bibetingelser, der er givet i Benth et al. [2007]. Den sørger for, at forwardprisen til tid 0 er lig spotprisen. Fælles for alle bibetingelserne er, at de er lineære i koe cienterne indeholdt i x. De kan derfor skrives på formen Ax = b. Når der som angivet i (8.8) integreres over udtrykket for
00 (t)2 kan kriteriefunktionen udtrykkes på matriceformen i ligning (8.15).
(8.15)
min x0 Hx x
Under bibetingelse af Ax = b og hvor matricen H er givet i ligning (8.16).
0 h1 .. H= . 0 hn
144 (5) ∆ 5 j
,
18∆(4) j (3) hj = 8∆j 0 0
(4)
(3)
0 0
18∆j
8∆j
(3) 12∆j
(2) 6∆j
0
(1)
(2)
6∆j
4∆j
0
0
0
0
0
0
0
(l)
∆j = tlj − tlj−1
0 0 0 0
(8.16)
(8.17)
Problemet omskrives til det ubegrænsede minimeringsproblem i (8.18) med Lagranges multiplikator metode.
min x0 Hx − λ0 (Ax − b) x,λ
56
(8.18)
Hvor løsningen [x∗ , λ∗ ] er løsningen til det lineære ligningssystem givet i (8.19).
0 2H A x 0 = A 0 λ b Således er koe cienterne i x∗ brugt til at beregne den udglattede kurve i gur 8.1.
57
(8.19)
8.2 Modi cerede simulationer Vi vil nu kalibrere simulationerne til den udglattede forwardkurve og derved danne modi cerede simulationer. Vores modi cerede simulationer er valgt som i (8.22). Her er det vigtigt at pointere, at de er valgt. De ikke er entydigt bestemt. Formen i (8.22) sikrer, at de rammer den udglattede forwardkurve i gennemsnit. Lad B være antallet af simulerede stier og Ptb være prisen for sti b ∈ {1, 2, ..., B} til tid t. b = P Pmod t B
Ptb
· f (t0 , t) b P /B b=1 t
B X
b Pmod t
= Pb=1 B
b=1
PB
PB
Ptb
b b=1 Pt
· f (t0 , t) · B
⇒
(8.20)
⇔
(8.21)
b Pmod t = f (t0 , t) B
(8.22)
b=1
50
Efter kalibrering
50
Før kalibrering
40
40
95% quantile 95% quantile
Mean
30
Mean
20
Eur/MWh
30
75% quantile
20
Eur/MWh
75% quantile
50% quantile
25% quantile
10
10
50% quantile
25% quantile
5% quantile
0
0
5% quantile
maj
jul
sep
nov
jan
mar
maj
jul
sep
nov
jan
mar
Figur 8.2: I det venstre panel ses de simulerede spot priser, før de er kalibreret til forwardkurven. I det højre ses de modi cerede simulationer, hvis gennemsnit til tid t er lig f (t0 , t).
58
Del II
Prisfastsættelse af en gaslagerkontrakt ved Least Squares Monte Carlo
59
9
Prisfastsættelse af gaslagerkontrakter
Dette afsnit er hovedsageligt baseret på artiklerne Boogert and De Jong [2008] og Boogert and de Jong [2011]. Først vil vi i afsnit 9.1 beskrive opbygningen af en gaslagerkontrakt. Dernæst angiver vi en række forsimplende antagelser i afsnit 9.2 og 9.3. Anvendelsen af Least Squares Monte Carlo metoden til værdiansættelsen gaslagere vil blive beskrevet i afsnittet 9.4 og vil blive illustreret i et simpelt eksempel i afsnit 9.5. Formålet med det illustrative eksempel er, at understøtte intuitionen bag LSMC metoden i prisfastsættelsen af gaslageret.
9.1 Gaslagerkontrakt Her betragtes en gaslagerkontrakt, der løber fra tid 0 til tid T + 1. Det er muligt at foretage handlinger på tidspunkterne t ∈ {1, 2, ..., T }. Tidsenheden måles her i dage. På en given dag kan man tilkøbe eller sælge gas. Derved pumpes der gas enten ud af eller ind i lageret. Man kan også undlade at handle, hvorved mængden af gas i lageret vil være uændret2 . Lad v(t) angive volumen af gas i lageret på starten af dag t og ∆v(t) være ændringen i volumen fra dag t til dag t + 1 dvs. ∆v(t) = v(t + 1) − v(t). Vi angiver volumen i MWh, der egentlig er en enhed for energi. Så når MWh benyttes skal det ses i lyset af, at der går 94.7265 m3 per MWh. Det antages endvidere, at prisen Pt , målt i EUR/MWh, er kendt ved begyndelse af dag t. Det antages endvidere, at priserne er likvide således, at handlerne, der foretages på dag t, ikke påvirker prisen. Volumen på dag t er givet ved summen af volumeændringerne fra tid 0 til tid t.
v(t) = v(0) +
t X i=1
2 Det
antages her, at volumen af gas ikke svinder.
60
∆v(i − 1)
(9.1)
Lad h(Pt , ∆v(t)) betegne payo ved en given handel på dag t. Dette er givet ved ligning (9.2).
h(Pt , ∆v(t)) =
−[(1 + a1 )Pt + b1 ]∆v(t) 0 −[(1 − a2 )Pt − b2 ]∆v(t)
∆v(t) > 0 injection ∆v(t) = 0 ingen handling
(9.2)
∆v(t) < 0 withdrawal
Her betragtes a1 og a2 som et mål for transaktionsomkostninger, mens b1 og b2 er et mål for bid-ask spread. Mængden af handlinger er begrænset af, at gaslageret har en maksimal og minimal kapacitetsgrænse for volumen af gas.
v min ≤ v(t) ≤ v max
(9.3)
Ligeledes er der en begrænsning på, hvor meget gas, der kan pumpes ud og ind af lageret på en enkelt dag. Den maksimale injection rate, dvs. den maksimale mængde, man kan pumpe ind, er givet ved imax . Man kan dog kun pumpe imax ind, hvis man ikke er begrænset af den maksimale kapacitet. Hvis spændet fra det nuværende niveau i volumen op til det maksimale volumen niveau, v max − v(t), er mindre end imax , er den maksimale injection v max − v(t). Tilsvarende betingelse gælder for withdrawal ved den nedre kapacitetsgrænse.
− min{wmax , v(t) − v min } ≤ ∆v(t) ≤ min{imax , v max − v(t)}
(9.4)
Værdien af en gaslagerkontrakt er nutidsværdien af de fremtidige betalingsstrømme, hvor der til hvert tidspunkt er foretaget et optimalt valg. Den optimale strategi π afhænger af tiden, gasprisen og volumen niveauet i lageret. π = {π(1, P1 (I1 ), v(1), I1 ), · · · , π(T, PT (IT ), v(T ), IT )}, hvor It angiver de faktorer, der har en indvirkning på gasprisen, og som er observerbare lige før beslutningen tages på tidspunkt t. Værdien ved til tid T +1 at have gas tilbage i lageret er givet ved skrotværdien q(PT +1 (IT +1 ), v(T + 1)). Det ønskes at maksimere værdien af lageret
61
ved at nde den optimale strategi. Dette giver anledning til følgende optimeringsproblem.
" sup EQ ∆v(t)
T X
! e−δt h(Pt (It ), ∆v(t))
# + e−δ(T +1) q(PT +1 (IT +1 ), v(T + 1))
(9.5)
t=1
Den optimale strategi må hverken bryde med begrænsningerne på kapaciteten i lageret eller begrænsningerne på ændringen i volumen niveau. Mængden af tilladte volumen niveauer på tid t de neres i ligning (9.6).
V(t) = {v | v min ≤ v ≤ v max }
(9.6)
Volumen niveauet i lageret på tid t indvirker på begrænsningen for mængden af tilladte ændringer. Ændringen i lageret må ikke medføre, at volumen niveauet falder under v min eller stiger over v max . Ændringen må ej heller overstige de maksimale injection og withdrawal rater. Mængden af tilladte handlinger angives ved D(t, v(t)).
D(t, v(t)) = {∆v | v min ≤ v(t) + ∆v ≤ v max ,
(9.7)
− wmax (t, v(t)) ≤ ∆v ≤ imax (t, v(t))}
9.2 Forsimplende antagelser En lang række markedsfriktioner og fysiske egenskaber ved operation af lageret vil blive underlagt en række forsimplende antagelser, da de ikke vurderes at spille en hovedrolle i prisfastsættelsen. Man kan fjerne disse antagelser, hvilket vil øge kompleksiteten af problemet.
• wmax (t, v(t)) og imax (t, v(t)) antages konstante. De antages således ikke at afhænge af volumen eller tid.
wmax (t, v(t)) = wmax imax (t, v(t)) = imax
62
• Det antages, at den gas, der ikke bliver solgt indenfor kontraktperioden på et år, er tabt. Det kan således ikke betale sig at have overskydende gas ved kontraktudløb. Det medfører, at man vil sælge gassen3 og ende med en volumen på v min . I ligning (9.5) og (9.16) vil det medføre q(PT +1 , v(T + 1)) = 0.
• Det antages, at v min = 0. • Transaktionsomkostningerne og bid-ask spreadet sættes til 0. Dvs. at a1 = a2 = b1 = b2 = 0. Payo-funktionen kan derfor reduceres til ligning (9.8).
h(Pt , ∆v(t)) = −Pt · ∆v(t)
(9.8)
9.3 Diskretisering af handlingsrum, tid og volume I dette afsnit vil der foretages en diskretisering af tiden, volumen og handlingsrummet. De tre variable kan i teorien antage uendeligt mange værdier, men reelt vurderes kun et endeligt antal muligheder. Mængden af uendelige valg skal derfor diskretiseres til en række endelige valg. Inddelingen af tid sker således, at der maksimalt kan træes et valg per dag. Handlingsrummet bliver diskretiseret, så der maksimalt bliver betragtet tre valg for et givent tidspunkt og med en given volumen. Før den nærmere speci cering af disse valg, vil vi beskrive en række specielle punkter (t, v(t)) i forhold til diskretiseringen i volumen niveauer. Kombinationen af tidspunkt og volumen niveau gør, at man kan være tvunget til at pumpe gas ud af lageret, hvis man skal kunne nå at tømme lageret. Givet antagelsen om, at tilbageværende gas i lageret ikke giver noget payo til tid T + 1, er det ønskværdigt at sælge alt gassen i lageret. Derfor er punkterne i (t, v(t)), hvor man lige akkurat kan nå at få tømt lageret til tid T + 1 specielle. De er specielle i den henseende, at man her reelt kun har et valg, nemlig at pumpe den maksimale mængde gas ud af lageret. Mængden af disse punkter benævnes 3 Det
vil give en positiv gevinst at sælge gassen, så længe prisen på gas overstiger transaktionsomkost-
ningerne.
63
som withdrawal-randen og er givet ved ligning (9.9).
wrand (t) = wmax · (T + 1 − t)
(9.9)
Det antages, at det ikke er ønskværdigt at have et volumen niveau, der overstiger wrand (t). Denne betingelse kombineres med begrænsningen v(t) ≤ v max . Det strengeste af de to krav, vil blive betegnet v rand , der er den øvre rand og er givet ved ligning (9.10). Den nedre rand er stadig givet ved v min .
v rand (t) = min{v max , wmax · (T + 1 − t)}
(9.10)
Her benyttes withdrawal-randen til at danne det diskrete grid i volumen niveauer. Der dannes således Nvolume = ceiling( v
max ·precision
wmax
) + 1 volumen matricer, hvor matricen indekseret
ved k ∈ {1, 2, ..., Nvolume − 1} har et volumen niveau på v(t, k) = (k · wmax )/precision og v(t, k = Nvolume ) = v max . I gur 9.1 er et eksempel med precision = 2 angivet, hvor withdrawal-randen er markeret. Volumen gridpunkterne dannes i algoritmen som Nvolume −1 punkter med ækvidistante afstande på wmax /precision samt det sidste punkt v(Nvolume ), der har en afstand, som er mindre end eller lig wmax /precision til det nærmeste underliggende nabopunkt. Dette er valget for diskretiseringen i tid og volumen.
64
v max
wmax
v min t1
t2
t3
t4
t5
t6
Figur 9.1: Grid i volumen og tid med withdrawal-randen markeret. Handlingsrummet diskretiseres ligeledes. Dette sker på en sådan måde, at der i et punkt
(t, v(t)) maksimalt vil vurderes tre valg, der som udgangspunkt er maksimal withdrawal, maksimal injection samt handlingen, hvor volumen niveauet ikke ændres. Der betragtes 5 typer af punkter i (t, v(t)), hvor forskellige handlinger tillades. 1. I punkterne
(ruder) i gur 9.2 er der kun et valg. Her sælges så meget som muligt,
og der udtages således den maksimale mængde fra lageret. Dette gør sig gældende for volumen niveauer, hvor man a) kun lige akkurat kan nå, eller b) ikke kan nå at pumpe gassen ud inden kontrakten udløber. Dette er tilfældet for (T + 1 − t) · wmax ≤ v(t) ≤
v max , hvor handlingen er givet ved ligning (9.11): ∆v(t) = − min{wmax , v(t) − v min }
2. I punkterne
(9.11)
• (cirkel) er det muligt at pumpe ind, ud eller gøre intet. Det ses i gur
9.2, at det i de øverste cirkel-punkter ikke vil være muligt at pumpe imax ind, men kun v max − v(t), da man rammer lagerets kapacitetsbegrænsning. Ligeledes vil det i de nederste cirkel-punkter ikke være muligt at pumpe wmax ud, men kun v(t) − v min . 65
Handlingerne som vurderes i cirkel punkterne er givet ved ligning (9.12).
min{imax , v rand (t + 1) − v(t)} ∆v(t) = 0 − min{wmax , v(t) − v min }
maksimal injection ingen handling
(9.12)
maksimal withdrawal
v max
wmax
v min t1
t2
t3
t4
t5
t6
Figur 9.2: Grid i volumen og tid med withdrawal-randen markeret.
3. I punkterne ( rkant) er handlingen begrænset til, at man kan gøre intet eller pumpe ud af lageret. Det er ikke muligt at pumpe ind, da volumen, v(t+1) = v(t)+∆v(t), enten vil bryde med kapacitetsbegrænsningen eller nå et niveau over withdrawal-randen. Derved er karakteriseret ved v(t) = v rand (t + 1) = min{wmax · (T − t), v max } med
t < T . Her vurderes derfor kun to handlinger, som er givet i ligning (9.13).
∆v(t) =
0
ingen handling
(9.13)
− min{wmax , v(t) − v min } maksimal withdrawal 4. I punkterne
N (trekant) er man tvunget til at pumpe ud, da man ved ikke at foretage 66
sig noget, ender over withdrawal-randen. Man er dog ikke tvunget til at pumpe den maksimale mængde ud af lageret. Den minimale mængde, der skal pumpes ud af lageret, er mængden, hvor volumen niveauet v(t + 1) ligger på withdrawal-randen. De to yderpunkter holdes op mod hinanden. Disse punkter er karakteriseret ved
v rand (t + 1) < v(t) ≤ v rand (t). Handlingerne, der vurderes, er givet ved ligning (9.14).
∆v(t) =
−(v(t) − v rand (t + 1))
minimal withdrawal
−wmax
maksimal withdrawal
(9.14)
5. I punkterne F (stjerne) kan man enten pumpe gas ind eller lade gasvolumen forblive den samme. Dette gør sig gældende, når v(t) = v min med t < T . Handlingen er givet ved ligning (9.15).
∆v(t) =
min{imax , v rand (t + 1) − v(t)}
maksimal injection
0
ingen handling
67
(9.15)
9.4 Least Squares Monte Carlo Værdien af en gaslagerkontrakt, der løber fra tid t til T + 1 angives ved U (t, Pt (It ), v(t)). Dette er værdien af kontrakten før handlingen ∆v(t), der foretages på dag t. På dag T + 1, hvor man ikke har ere handlinger, er værdien af kontrakten givet ved skrotværdien af den tilbageværende mængde gas i lageret.
U (T + 1, PT +1 (IT +1 ), v(T + 1)) = q(PT +1 (IT +1 ), v(T + 1))
(9.16)
På dag T foretages den sidste handling. Her er værdien af lagerkontrakten givet ved cash owet for den sidste handling plus den tilbagediskonterede værdi af kontrakten på tid
T + 1.
U (T, PT (IT ), v(T )) =
max ∆v∈D(T,v(T ))
EQ h(PT (IT ), ∆v(T )) + e−δ U (T + 1, PT +1 (IT +1 ), v(T + 1)) (9.17)
Dette er starten på en iterativ løsningsmetode, da værdien af kontrakten til et generelt tidspunkt t derved dels består af cash owet h(Pt (It ), ∆v(t)), der knytter sig til handlingen på den givne dag, samt nutidsværdien af handlingerne i den resterende kontraktperiode. Denne fortsættelsesværdi betegnes ved C(t, Pt (It ), v(t), ∆v(t)) og vil refereres til som continuation value. Det ses, at continuation value afhænger både af volumen niveauet og ændringen i volumen. Boogert and De Jong [2008] påpeger, at continuation value på tidspunkt t kun afhænger af det volumen niveau, der opnås efter ændringen ∆v(t), dvs. niveauet til næste tidspunkt v(t + 1) = v(t) + ∆v(t). Continuation value kan så skrives som C(t, Pt (It ), v(t) +
∆v(t)) = C(t, Pt (It ), v(t+1)). Continuation value de neres ved relationen i ligningen (9.18). C(t, Pt (It ), v(t + 1)) = EQ e−δ U (t + 1, Pt+1 (It+1 ), v(t) + ∆v(t))
(9.18)
For et generelt t i kontraktperioden er værdien af kontrakten givet ved det dynamiske optimeringsproblem i ligningerne (9.16) og (9.19). Problemet løses rekursivt ved baglæns
68
induktion.
U (t, Pt (It ), v(t)) =
max
{h(Pt (It ), ∆v(t)) + C(t + 1, Pt+1 (It+1 ), v(t + 1))}
∆v∈D(t,v(t))
(9.19)
Det ønskes at nde en optimal beslutningsregel i det dynamiske optimeringsproblem. For at træe optimale valg skal et estimat for continuation value ndes. Dette kan gøres ved brug af Least Squares Monte Carlo. På tid t approksimeres continuation value ved en lineær sum af q ∈ {1, ..., Q} basisfunktioner φq (·) af en række tilstandsvariable X .
C(t, X) ≈
Q X
φq (t, X)βq
(9.20)
q=1
Tilstandsvariablene kunne som i (9.18) være X = [ t
Pt (It ) v(t + 1) It ]. Der estimeres
et sæt af β -koe cienter for hvert punkt i griddet (t, v(t + 1)) beskrevet i afsnit 9.3. Derfor undlades det at inkludere basisfunktioner med tilstandvariablene v(t + 1) og t. Det antages, at continuation value for hver af de B uafhængige simulerede stier approksimativt er givet ved ligningen (9.21). b C b (t, Ptb (It ), v(t + 1)) ≈ e−δ Y b (t + 1, Pt+1 (It+1 ), v(t + 1))
(9.21)
b Her betegner Y b (t + 1, Pt+1 (It+1 ), v(t + 1)) den akkumulerede fremtidige værdi af de penge-
strømme, der opnås for simulation b, hvor
• Den optimale beslutningsregel følges fra tid t + 1 til tid T . • Der på tidspunkt t + 1 er et volumen niveau på v(t + 1) og prisstien er simuleret til b Pt+1 (It+1 ).
Estimaterne βˆq estimeres ved mindste kvadraters metode ud fra estimationsligningen (9.22).
b e−δ Y b (t + 1, Pt+1 (It+1 ), v(t + 1)) =
Q X q=1
69
φq (Pt (It ))βq (t, v(t + 1)) + bt
(9.22)
Hvor βˆq (t, v(t + 1)) her implicit angiver beslutningsreglen. Estimatet for continuation value ndes som:
ˆ P b (It ), v(t + 1)) = C(t, t
Q X
φq (Ptb (It ))βˆq (t, v(t + 1))
(9.23)
q=1
Hvis en handling medfører et skift til et volumen niveau, der ikke er indeholdt i volumen griddet, da ndes continuation value for dette volumen niveau ved en lineær interpolation af continuation value mellem de to nabopunkter.
9.4.1
Basisfunktioner
I LSMC approksimeres til hvert tidspunkt en continuation value, der benyttes til at lave den optimale beslutningsregel. Denne approksimation baseres på en lineær sum af basisfunktioner, som vist i ligning (9.20). Mange forskellige typer af polynomier kan benyttes som basisfunktion. I denne opgave vil det undersøges, hvilken eekt valget af polynomium har på prisen af lageret. Mere speci kt vil re typer af polynomier blive afprøvet. En simpel sum af monomials, Laguerre, Weighted Laguerre og Legendre polynomier. Disse re polynomier kan opskrives rekursivt og er vist i tabel 9.1.
1. medlem
M0 (x) = 1
3. medlem
−x
W L0 (x) = e 2
L2 (x) =
−x
W L1 (x) = e 2
k 'te medlem
M2 (x) = x2
L1 (x) = 1 − x
L0 (x) = 1
LE0 (x) = 1
2. medlem
M1 (x) = x
· (1 − x)
LE1 (x) = 2x − 1
W L2 (x)
x2 −4x+2 2 −x 2 = e 2 · ( x −4x+2 ) 2 2
LE2 (x) = 6x − 6x + 1
Mk+1 (x) = xk+1
k ·L Lk+1 (x) = 2k+1−x · Lk (x) − k+1 k−1 (x) k+1 −x
W Lk+1 (x) = e 2 LEk+1 (x) =
2k+1 k+1
k ·L · ( 2k+1−x · Lk (x) − k+1 k−1 (x)) k+1 k ·P · xPk (x) − k+1 k−1 (x)
Tabel 9.1: Mk (x) angiver monomials med medlem k for statevariablen x, Lk (x) angiver tilsvarende Laguerre, W Lk (x) Weigthed Laguerre og LEk (x) Legendre. For de tre sidstnævnte polynomier gælder, at de er ortogonale over forskellige intervaller. Denne egenskab kan ifølge Areal et al. [2008] være en fordel, når der skal estimeres koe cienter til continuation value, da man i højere grad undgår multikollinearitet i designmatricen. Man afhjælper herved situationen, hvor designmatrice bliver tilnærmelsesvis singulær, og koe cienter enten bliver bestemt med stor usikkerhed eller ikke kan beregnes. I testkørsler med få antal simulationer har det i udpræget grad betydning, hvilken type af polynomium, der benyttes. Ved et større antal simulationer bør prisen være den samme for Laguerre, Legendre og Monomials. Dette skyldes, at de tre polynomier kan opskrives som en linear 70
kombinationer af hinanden. Weighted Laguerre falder udenfor denne kategori grundet vægten. Dette resulterer i, at vi ikke forventer den samme pris ved brug af dette polynomium, som ved brug af de andre. I gur 9.3 ses syv medlemmer af Monomials. I dette og de næste plots er det at foretrække, hvis medlemmerne ikke er divergeret i det interval, hvor priserne hovedsageligt be nder sig. Dette vil sikre velbestemte estimater. Priserne be nder sig i intervallet fra 5 til 70, og
6e+05 4e+05
M1(x) M2(x) M3(x) M4(x) M5(x) M6(x) M7(x)
0e+00
2e+05
Funktionsværdi
8e+05
1e+06
monomials er derfor måske ikke det optimale polynomium ved et lavt antal simulationer.
0
20
40
60
80
100
Gaspris
Figur 9.3: Monomials. Weighted Laguerre tilhører de ortogonale polynomier. Denne er vist i gur 9.4 og er ortogonal i intervallet [0; ∞[. På guren ses det, at Weighted Laguerre ikke er divergeret i prisintervallet fra [0; 70], men til gengæld konvergerer mod nul. Dette sker ved en pris på omkring 40, hvilket kan medføre, at Weighted Laguerre performer dårligere for priser over 40. En mulig løsning er at skalere sine priser, så priserne ligger i intervallet [0; 40]. Dette er undladt.
71
1.0 0.5 0.0
Funktionsværdi
−1.0
−0.5
1.0 0.5 0.0
Funktionsværdi
−0.5 −1.0
WL1(x) WL2(x) WL3(x) WL4(x) WL5(x) WL6(x) WL7(x)
0
5
10
15
20
25
0
20
40
Gaspris
60
80
100
Gaspris
Figur 9.4: Weighted Laguerre. I gur 9.5 er syv medlemmer af Laguerre plottet. Grafen til venstre viser, hvordan Laguerre medlemmerne bevæger sig for priser under 25. På grafen til højre er tilsvarende plottet for priser under 100. Her ses det, at medlemmerne, som er større end L3 (x), relativ hurtigt divergerer. For at sikre velbestemte estimater vil der i analysedelen benyttes et stort antal simulationer. Det samme gælder for Legendre, der ses plottet i gur 9.6. For priser under 100 divergerer medlemmerne af Legendre. I grafen til venstre er Legendre plottet i intervallet [−1, 1]. Det er her vigtigt at pointere, at gasprisen selvfølgelig ikke kan være negativ. Grafen viser intervallet, hvor Legendre er ortogonal. For at få det maksimale ud af at benytte Legendre skal priserne skalleres ned, så de be nder sig i intervallet [−1, 1]. Dette er undladt, da Areal et al. [2008] som sagt viser, at prisen bliver velbestemt ved et stort antal
0
5
10
15
20
25
5e+05
L1(x) L2(x) L3(x) L4(x) L5(x) L6(x) L7(x)
0e+00 −1e+06
Funktionsværdi
100 0 −100 −300
Funktionsværdi
200
300
1e+06
simulationer.
0
Gaspris
20
40
60 Gaspris
Figur 9.5: Laguerre.
72
80
100
1e+06 −1e+06
0.0
0.5
5e+05
Funktionsværdi −0.5
0e+00
1.0 0.5 0.0
Funktionsværdi
−0.5 −1.0 −1.0
LE1(x) LE2(x) LE3(x) LE4(x) LE5(x) LE6(x) LE7(x)
1.0
0
20
Gaspris
40
60
80
100
Gaspris
Figur 9.6: Legendre.
9.5 Illustrativt eksempel I dette afsnit vil LSMC metoden gennemgås ved brug af et illustrativt eksempel. I eksemplet vil vi speci kt vise, hvordan typen af valg afhænger af kombinationen af volumen niveau og tidspunkt. Parameterværdierne er sat således, at alle scenarierne fra gur 9.2 dækkes. Følgende inputs er valgt: i = 1.65 MWh/dag, w = 2.20 MWh/dag, d = 0.9967, V max =
3 MWh, V min = 0 MWh og precision = 2. I gur 9.7 ses griddet i tid og volumen. 3.0 MWh 2.2 MWh 1.1 MWh 0.0 MWh
t1
t2
t3
t4
Figur 9.7: Gridpunkter for det illustrative eksempel. I det følgende afsnit vil der både refereres til de procentvise og nominelle volumen niveauer. Det betyder, at v(t, k = 4) = 3 MWh svarer til et volumen niveau på 100%. Ligeledes vil
v(t, k = 3) = 2.2 MWh svare til 73.3%, v(t, k = 2) = 1.1 MWh til 36.7% og v(t, k = 1) = 0 MWh svarer til et volumen niveau på 0%. I tabel 9.2 ses de 8 simulerede prisstier. Alle de simulerede stier starter ved Pt0 = 21.4 EUR/MWh. 73
Simulerede t1 t2 21.0 15.8 21.0 17.2 21.1 12.3 19.8 14.7 20.7 15.4 20.0 24.7 19.5 16.4 19.8 13.5
prisstier t3 t4 8.9 8.5 22.6 21.1 10.7 7.8 18.4 19.0 21.0 27.5 22.1 19.9 18.2 31.0 19.1 26.8
Tabel 9.2 Først fyldes de re payo matricer op på sluttidspunktet t4 , hvor der kun er en reel handling, hvilket er at pumpe den maksimale mængde gas ud af lageret. Dette ses i gur 9.7 ved, at alle punkterne til tid t4 er ruder. Der pumpes således en mængde ud på ∆v = −min(wmax , v(t)−
vmin ), hvilket svarer til 2.2 MWh for både v(t4 , k = 3) og v(t4 , k = 4). For v(t4 , k = 1) og v(t4 , k = 2) er det ikke muligt. Der pumpes derfor henholdsvis intet og 1.1 MWh ud af lageret. Det giver payo matricerne i tabel 9.3. Volumen = 0.0 MWh t1 t2 t3 t4 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
Payo matricer Volumen = 1.1 MWh Volumen = 2.2 MWh t1 t2 t3 t4 t1 t2 t3 t4 9.4 18.8 23.2 46.4 8.6 17.1 20.9 41.9 30.3 60.6 21.8 43.7 34.2 68.3 29.5 59.0
Volumen = 3.0 MWh t1 t2 t3 t4 18.8 46.4 17.1 41.9 60.6 43.7 68.3 59.0
Tabel 9.3 Næste skridt er at bestemme continuation values. Disse ndes ved at benytte ligning (9.21). Først ndes de tilbagediskonterede fremtidige payos. Til dette benyttes d, og resultatet ses i tabel 9.4.
74
Tilbagediskonterede payos t3 0.0 MWh 1.1 MWh 2.2 MWh 3.0 MWh 0.0 9.3 18.7 18.7 0.0 23.1 46.2 46.2 0.0 8.5 17.1 17.1 0.0 20.9 41.8 41.8 0.0 30.2 60.4 60.4 0.0 21.8 43.5 43.5 0.0 34.0 68.1 68.1 0.0 29.4 58.8 58.8 Tabel 9.4 Estimatet for den fremtidige værdi af gaslagerkontrakten ndes ved at benytte lineær regression. Her har vi valgt monomials med gasprisen som tilstandsvariabel. Denne indgår til og med anden orden, hvilket svarer til et andengradspolynomium af priserne. Disse regresseres på de tilbagediskonterede payos angivet ved Y .
eδ · Y (t + 1, Pt+1 , v(t + 1; k)) ≈ β0 (t, k) + β1 (t, k) · Pt + β2 (t, k) · Pt2
(9.24)
Regressionskoe cienter er angivet i tabel 9.5.
βˆ0 βˆ1 βˆ2
Regressionskoe cienter t3 0.0 MWh 1.1 MWh 2.2 MWh 0.00 -53.63 -107.27 0.00 8.79 17.58 0.00 -0.24 -0.47
3.0 MWh -107.27 17.58 -0.47
Tabel 9.5 Koe cienterne skal benyttes til at estimere continuation value for hver enkelt sti. Der er i alt B stier, hvor b er indekset for en speci k sti. b ∈ {1, ..., B}.
Cˆ b (t, Ptb , v(t + 1, k)) = β0 (t, k) + β1 (t, k) · Ptb + β2 (t, k) · (Ptb )2
(9.25)
I ligning (9.26) er continuation value udregnet for den første simulation med et volumen
75
niveau på 1.1 MWh.
Cˆ 1 (t3 , 8.9304, 1.1) = −53.63 + 8.7899 · 8.9304 − 0.2373 · 8.93042 = 5.9422
(9.26)
På samme måde ndes de øvrige continuation values for tidspunktet t3 . Disse er angivet i tabel 9.6. Tallet markeret med orange er den udregnede continuation value fra ligning (9.26). Continuation value t3 0.0 MWh 1.1 MWh 2.2 MWh 3.0 MWh 0.0 5.9 11.9 11.9 0.0 23.7 47.5 47.5 0.0 13.2 26.5 26.5 0.0 27.8 55.5 55.5 0.0 26.3 52.7 52.7 0.0 24.8 49.5 49.5 0.0 27.7 55.5 55.5 0.0 27.7 55.4 55.4 Tabel 9.6 I hvert (t, k)-punkt ønskes det at træe det optimale valg, der maksimerer ligning (9.19). I ligning (9.27) er værdien af kontrakten speci ceret til tid t for de tre diskrete valg, injection, gøre ingenting og withdrawal.
Uinj (t, k)b = −Ptb · vinj + Cˆ b (t, Ptb , v(t, k) + vinj ) Unull (t, k)b = Cˆ b (t, Ptb , v(t, k))
(9.27)
Uwdraw (t, k)b = Ptb · vwdraw + Cˆ b (t, Ptb , v(t, k) − vwdraw ) For at eksempli cere, hvordan det optimale valg træes, udregnes nedenfor værdien af kontrakten for sti nummer 1 ved alle re volumen niveauer. Der ses på de re punkter til tidspunkt t3 , som er markeret med orange i gur 9.8.
76
3.0 2.2 1.1 0.0
t1
t2
t3
t4
Figur 9.8: Gridpunkter for det illustrative eksempel. Først undersøges punktet F. I dette punkt er withdrawal ikke en mulighed, da der ikke er noget gas i lageret. Dvs. at det optimale valg enten er at pumpe gas ind i lageret eller gøre intet. Hvis der pumpes gas ind i lageret, opnås et volumen niveau på 1.65 MWh. Dette ligger mellem grid volumen niveauerne 1.1 MWh og 2.2 MWh, der er indeholdt i griddet. Continuation value til næste tidspunkt af at stå med en volumen på 1.65 MWh, ndes derfor ved lineær interpolation mellem continuation values i gridpunkterne 1.1 MWh og 2.2 MWh. I den lineære interpolation angiver xi vægten på det øvre volumen niveau ved injection, dvs. her vægten på continuation value ved en volumen på 2.2 MWh.
(v[t3 , 1] + min(i, v max − v[t3 , 1]) − v[t3 , 2]) xi = v[t3 , 3] − v[t3 , 2] (0 + 1.65 − 1.1) = = 0.5 2.2 − 1.1 Continuation value for et volumen niveau på 1.65 MWh udregnes ved at vægte continuation values fra tabel 9.6.
Cˆ 1 (t3 , 8.9304, 0 + 1.65) = 0.5 · 11.88348 + (1 − 0.5) · 5.941738 = 8.912609
(9.28)
Ligning (9.27) kan nu opstilles speci kt for F. Det huskes fra tabel 9.6, at continuation
77
value til tid t3 med et tomt lager, selvfølgelig er 0.
Uinj (t3 , k = 1) = −8.9304 · 1.65 + 8.912609 = −5.822551
(9.29)
Unull (t3 , k = 1) = 0 Det optimale valg i dette punkt er derfor ikke at pumpe gas ind i lageret. Næste skridt er nde det optimale valg i punktet •. Volumen niveauet i punktet er 1.1 MWh, hvilket er mindre end den maksimale withdrawal rate på 2.2 MWh. Derved kan der maksimalt pumpes 1.1 MWh ud af lageret. Pumpes en mængde på 1.65 MWh ind i lageret resulterer det i en mængde på 2.75 MWh. Det ønskes ikke at fylde lageret, så volumen niveauet overstiger v rand fra (9.10), hvilket er tilfældet, hvis der injectes 1.65 MWh. Til næste tidspunkt kan man kun nå at pumpe 2.2 MWh ud af lageret, og man bliver derfor ikke belønnet for de resterende 0.55 MWh. Derved ønskes det maksimalt at pumpe 1.1 MWh ind i lageret. Dette medfører, at det ikke er nødvendigt med lineær interpolation. Både ved withdrawal og injection rammes punkterne i griddet. Ligningerne i (9.30) viser værdien af kontrakten ved de tre valg.
Uinj (t3 , k = 2) = −8.930425 · 1.1 + 11.88348 = 2.06004 Unull (t3 , k = 2) = 5.941738
(9.30)
Uwdraw (t3 , k = 2) = 8.930425 · 1.1 + 0 = 9.82344 Den højeste værdi af kontrakten opnås derfor i dette punkt ved withdrawal. Punktet er specielt på den måde, at det ikke kan betale sig at pumpe gas ind i lageret. Ved injection vil v rand blive brudt. De to mulige valg er derfor at gøre intet eller at pumpe gas ud af lageret. Volumen i lageret er 2.2 MWh, hvilket er lig den maksimale withdrawal. I dette tilfælde vil man derfor pumpe 2.2 MWh ud af lageret, hvis det optimerer værdien af kontrakten.
Unull (t3 , k = 3) = 11.88348 Uwdraw (t3 , k = 3) = 8.930425 · 2.2 + 0 = 19.64688 I kan det derfor bedst betale sig at pumpe 2.2 MWh ud af lageret. 78
(9.31)
Det sidste orange punkt er
N. I dette punkt ønskes det at pumpe gas ud af lageret, da man
både ved at pumpe gas ind i lageret, eller ved at gøre ingenting opnår et volumen niveau over v rand . Det resulterer i de to mulige valg, hvor der henholdsvis pumpes lidt eller meget ud af lageret. Volumen i punktet er 3 MWh og det mindste, der skal pumpes ud for ikke at bryde v rand er 0.8 MWh. Den anden mulighed er at pumpe wmax = 2.2 MWh ud af lageret. I dette tilfælde ender volumen niveauet på 0.8 MWh. Dette volumen niveau er ikke indeholdt i griddet, hvorfor der interpoleres mellem de nærliggende volumen niveauer på 0 MWh og
1.1 MWh. I den lineære interpolation er xw vægten på den øvre volumen ved withdrawal, dvs. her er den vægten på continuation value ved en volumen på 1.1 MWh.
(v[t3 , 4] − min(w, v[t3 , 4] − v min ) − v[t3 , 1]) xw = v[t3 , 2] − v[t3 , 1] (3 − 2.2 − 0) = = 0.72727 2.2 − 1.1 Continuation value for et volumen niveau på 0.8 MWh kan ndes ved at vægte continuation values fra tabel 9.6.
Cˆ 1 (t3 , 8.9304, 3 − 2.2) = 0.72727 · 5.941738 + (1 − 0.72727) · 0 = 4.321264 Ligningerne i (9.32) viser værdien af kontrakten for hvert af de to valg.
Uwdrawmin (t3 , k = 4) = 8.930425 · 0.8 + 11.88348 = 19.027800
(9.32)
Uwdrawmax (t3 , k = 4) = 8.930425 · 2.2 + 4.321264 = 23.968144 Endnu en gang er resultatet, at der skal pumpes den maksimalt mulige mængde ud af lageret. I tabel 9.7 er værdien af kontrakten udregnet for samtlige volumen niveauer, stier og handlinger.
79
Volumen = 0.0 MWh w 0 i | 0.0 -5.8 0.0 -1.7 0.0 2.2 0.0 11.3 0.0 4.9 0.0 0.7 0.0 11.6 0.0 10.1
Volumen = w 0 9.8 5.9 24.9 23.7 11.8 13.2 20.3 27.8 23.1 26.3 24.3 24.8 20.0 27.7 21.0 27.7
Beslutningsmatrice 1.1 MWh Volumen = 2.2 MWh i w 0 i 2.1 19.6 11.9 | 22.6 49.8 47.5 14.7 23.5 26.5 35.3 40.5 55.5 29.6 46.2 52.7 25.2 48.6 49.5 35.5 40.0 55.5 34.4 42.0 55.4
Volumen = 3.0 MWh wmax wmin i 24.0 19.0 | 67.1 65.6 33.2 35.1 60.7 70.3 65.3 69.5 66.6 67.2 60.2 70.1 62.1 70.7
Tabel 9.7 Værdierne i tabel 9.7 udgør beslutningsgrundlaget, hvor valget, som har den største værdi, træes. I tabel 9.8 er de optimale valg tabelleret for alle stier og volumen niveauer.
0.0 MWh 0 0 i i i i i i
Beslutning til tid t3 1.1 MWh 2.2 MWh w w w w i 0 i 0 i 0 i 0 i 0 i 0
3.0 MWh wmax wmax wmin wmin wmin wmin wmin wmin
Tabel 9.8 Næste step er at opdatere payo matricerne til tid t4 . Dette gøres baseret på valgene, der træes i 9.8. Som eksempel kan ses på sti nr. 3 i situationen, hvor volumen niveauet er nul på tid t3 . Her er det vurderet optimalt at pumpe ind i lageret. Dette gør at volumen niveauet for sti 3 på tid t4 bliver 1.65 MWh. Værdien på tid t4 af at have 1.65 MWh i lageret er så ikke længere nul som vist i tabel 9.3, men i stedet 1.65 MWh gange prisen på
7.8 EUR/MWh. Et andet eksempel på, hvorledes matricen opdateres, er givet for den første sti ved et volumen niveau på 3. Det optimale valg er her at pumpe den maksimale mængde ud af lageret. Payomatricerne til tid t3 og t4 opdateres på følgende vis:
80
Payot3 ,k=4 = wmax · Pt13 = 2.2 · 8.9304 = 19.6469 Payot4 ,k=4 = (v(t3 , k = 4) − wmax ) · Pt14 = (3 − 2.2) · 8.5245 = 6.8196 Først opdateres payo matricen til tid t3 , hvor det vides fra tabel 9.8, at den optimale handling er at pumpe wmax ud af lageret. Det giver et payo på 19.6469 EUR. Den optimale handling til tid t4 er at pumpe så meget ud af lageret som muligt. Den resterende mænge gas i lageret er 0.8 MWh, som ganges på prisen til t4 . Det giver payoet på de 6.8196 EUR. De opdaterede payo matricer ses i tabel 9.9, hvor ovenstående payo-eksempel er markeret med orange. Volumen = 0.0 MWh t1 t2 t3 t4 0.0 0.0 0.0 0.0 -17.7 12.8 -30.4 31.4 -34.6 45.4 -36.4 32.8 -30.0 51.2 -31.5 44.2
Opdaterede payo matricer Volumen = 1.1 MWh Volumen = 2.2 MWh t1 t2 t3 t4 t1 t2 t3 t4 9.8 0.0 19.6 0.0 24.9 0.0 49.8 0.0 -11.8 17.1 0.0 17.1 -20.3 41.9 0.0 41.9 -23.1 60.6 0.0 60.6 -24.3 43.7 0.0 43.7 -20.0 68.3 0.0 68.3 -21.0 59.0 0.0 59.0
Volumen = t1 t2 t3 19.6 49.8 8.6 14.7 16.8 17.7 14.6 15.3
3.0 MWh t4 6.8 16.9 17.1 41.9 60.6 43.7 68.3 59.0
Tabel 9.9 I det ovenstående afsnit er der foretaget udregninger for t3 . Samtlige skridt skal ligeledes gennemgås for t2 og t1 . Her er nedenstående fremgangsmåde benyttet. 1. Udregn nutidsværdien af fremtidige payos. 2. Regresser nutidsværdien på tilstandsvariablene. 3. Find continuation values ved at tte regressionen. 4. Lav beslutningsmatricer ud fra de maksimalt 3 valg. 5. Benyt de fundne optimale valg til at opdatere payo matricerne. 81
Når ovenstående fremgang er foretaget for tid t2 og t1 er payomatricerne fyldte. Disse er vist i tabel 9.10.
t1 0.0 0.0 0.0 -32.6 0.0 -33.0 0.0 -32.7 t1 46.2 46.1 46.4 43.5 45.5 0.0 43.0 0.0
Volumen t2 -26.0 -28.4 -20.3 -18.0 -25.4 40.8 -27.0 -16.5 Volumen t2 -26.0 -28.4 -20.3 -24.3 -25.4 54.4 -27.0 -10.8
Opdaterede payo = 0.0 MWh t3 t4 14.7 0.0 37.4 0.0 -5.9 17.1 12.4 41.9 -11.5 60.6 -36.4 32.8 -10.0 68.3 12.9 59.0 = 2.2 MWh t3 t4 14.7 0.0 37.4 0.0 -5.9 17.1 -10.1 41.9 -11.5 60.6 -36.4 32.8 -10.0 68.3 15.3 59.0
matricer til tid t1 Volumen = 1.1 MWh t1 t2 t3 t4 23.1 -26.0 14.7 0.0 23.1 -28.4 37.4 0.0 23.2 -20.3 -5.9 17.1 0.0 -24.3 10.1 41.9 22.7 -25.4 -11.5 60.6 0.0 27.2 -36.4 32.8 21.5 -27.0 -10.0 68.3 0.0 -22.3 10.5 59.0 Volumen = 3.0 MWh t1 t2 t3 t4 46.2 -26.0 18.3 3.4 46.1 -28.4 46.4 8.4 46.4 -20.3 2.7 17.1 43.5 -24.3 4.6 41.9 45.5 -25.4 5.2 60.6 43.9 19.8 -36.4 32.8 43.0 -27.0 4.5 68.3 43.6 -22.3 4.8 59.0
Tabel 9.10 Volumen niveauerne i payo marticerne angiver nu start volumen i lageret til tid t1 , da der ikke kan foretages handlinger til tid t0 . Det vil sige, den første matrice med en volumen på
0.0 MWh, beskriver payo ved at starte med intet i lageret og slutte med et tomt lager. Ligeledes beskriver volumen niveauet på 3.0 MWh situationen, hvor der startes med et fyldt lager og sluttes med tomt lager. LSMC strategien sammenholdes nu med en simpel alternativ strategi. Strategien går i alt sin enkelthed ud på, at fylde lageret helt op for derefter at tømme det så sent som muligt. Det antages, at der startes med et tomt lager. Grundet setupet i eksempelet, betyder det, at der til t1 og t2 pumpes gas ind i lageret, mens der til t3 og t4 pumpes gas ud af lageret.
∆v(t) er speci ceret i tabel 9.11.
82
∆v(1) 1.65
∆v(2) 1.35
∆v(3) -0.8
∆v(4) -2.2
Tabel 9.11: ∆v(t) ved den alternative strategi. Payo udregnes for samtlige stier ved at benytte ligning (9.8). Payo matricen for den alternative strategi er vist i tabel 9.12. -34.67 -34.59 -34.79 -32.65 -34.11 -32.95 -32.22 -32.72
Volumen = 0.0 MWh -21.29 7.14 -23.24 18.12 -16.63 8.56 -19.86 14.74 -20.82 16.79 -33.39 17.67 -22.13 14.55 -18.21 15.27
18.75 46.35 17.11 41.90 60.58 43.68 68.31 58.97
Tabel 9.12 Det bemærkes, at for sti nummer re og otte, er den optimale handling den samme for både LSMC strategien og den alternative strategi. Her kan det bedst betale sig at pumpe gas ind i lageret til t1 og t2 , og pumpe det ud til t3 og t4 . Den opmærksomme læser vil ligge mærke til, at værdien af payo dog ikke er helt ens. Sti 4 8
-32.6 -32.7
LSMC -18.0 12.4 -16.5 12.9
41.9 59.0
-32.6 -32.7
Alternativ strategi -19.9 14.7 -18.2 15.3
41.9 59.0
Tabel 9.13 Dette skyldes, at gridpunkterne i dette simple eksempel ikke er nmaskede nok. Hvis en handling medfører, at et volumen niveau uden for griddet rammes, så interpoleres der mellem nabopunkterne. Speci kt påvirker det ovenstående payo i tabel 9.13 ved, at der til t1 injectes 1.65. Dette punkt ligger midt mellem 1.1 MWh og 2.2 MWh, hvilket resulterer i en vægt på 0.5 på nabopunkterne med volumen niveauerne 1.1 og 2.2. Til t2 er den optimale handling for disse to stier igen at injecte. Det medfører, at der injectes (3 MWh −
1.65 MWh) = 1.35 MWh for at fylde lageret helt op. I LSMC ndes den optimale injection 83
som sagt ved en vægtning af 1.1 MWh og 2.2 MWh. I punktet 1.1 MWh er det ikke muligt at fylde lageret helt op. Det er her kun muligt at fylde lageret til et niveau på 2.75 MWh. Når der interpoleres, så vægtes payo i punktet 1.1 MWh med 0.5, hvilket gør, at payo er beregnet ud fra den gennemsnitlige injection rate, som er 1.225 MWh i stedet for 1.35 MWh. (9.33)
0.5 · min(2.2 + 1.65, 3) + 0.5 · min(1.1 + 1.65, 3) − 1.65 = 1.225
Løsningen til ovenstående problematik er at gøre griddet tættere. Dette sker i praksis ved at vælge en højere precision. I eksemplet er denne sat til 2 og allerede ved en precision på 3 er griddet så nmasket, at eekten af interpolationen er forsvundet. I gur 9.9 er det akkumulerede payo for de 8 stier plottet for både LSMC og den alternative strategi. Figuren giver et indtryk af, hvordan der investeres i gaslageret de første perioder for at kunne skabe overskud i sidste ende. Det mest interessante ved dette plot er selvsagt værdien af lageret til tid t4 . I tabel 9.14 er payo for hver sti nutidsværdien til tid t0 ,
30 20 10 0 −10 −70
−50
−30
Samlet værdi af lagret
0 −10 −30 −70
−50
Samlet værdi af lagret
10
20
30
hvorefter de er summeret.
t0
t1
t2
t3
t4
t0
t1
t2
t3
t4
Figur 9.9: Det akkumulerede payo for de 8 stier. LSMC til venstre og den alternative strategi til højre. 84
Simulation 1 2 3 4 5 6 7 8 Gennemsnit
LSMC -11.25 8.78 -9.13 3.22 23.08 3.94 30.64 21.93 8.90
Alternativ strategi -30.12 6.11 -25.83 3.68 21.73 -5.41 27.72 22.62 2.56
Tabel 9.14 Gennemsnitsværdien af lageret ved de to strategier er vist nederst i tabel 9.14. Det forventede payo af at benytte LSMC er væsentligt større end ved den alternative strategi. Det skal dog huskes, at det er et meget lille illustrativt eksempel, der er lavet for at vise, hvordan LSMC fungerer i gaslager setuppet.
9.5.1
Rentens påvirkning på LSMC
I dette afsnit vil det illustrative eksempel benyttes til at vise, hvilken eekt en renteændring kan have på prissætningen af lageret i LSMC. Intuitivt vil en stigning i renten medfører, at fremtidige payos bliver tilbagediskonteres hårdere. Isoleret set vil dette mindske prisen på lageret, da man må købe gassen før, at man kan sælge den. Stiger renten kan beslutningsreglen ændre sig. Antag, at der indgås en handel, hvor der købes 1 MWh gas i dag med en omkostning på 24 EUR, der så sælges igen om 90 dage til en forventet pris på 24.10 EUR for den pågældende MWh. Så vil man ved en årlig rente på 1% have en positiv forventet nutidsværdi af den samlede handel på 0.04 EUR. Derimod vil man ved en rente på 2% have en negativ forventet nutidsværdi på 0.019 EUR. Således vil man ved de to forskellige renteniveauer træe forskellige beslutninger. Når det er sagt, så er diskonteringseekten meget lille specielt, når handlerne indgås dagligt. For at kunne handle optimalt, udregnes continuation values til næste tidspunkt for de pågældende handlinger. Hvis værdien af fremtidige payos mindskes grundet diskonteringen, vil man i visse tilfælde fravælge at pumpe gas ind i lageret eller sælge gassen i lageret hurtigere. 85
Vi vil nu vise, hvordan renten påvirker beslutningsregelen. Der tages udgangspunkt i et lager, der ved starttidspunktet har et volumen niveau på 0 %. Det illustrative eksempel er genkørt med renter gående fra 0 % til 5 % i hop af 1 procentpoint. Der sker en ændring i handelsstrategien for sti nummer 4 ved springet i renten fra 1 % til 2 %. Payo matricerne for de to rentesatser er vist i tabel 9.15. Opdaterede payo matricer til tid t1
t1 0.00 0.00 0.00 -32.65 0.00 -32.95 0.00 -32.72
Rente t2 -26.02 -28.40 -20.33 -18.02 -25.45 40.81 -27.05 -16.52
= 1% t3 14.74 37.36 -5.88 12.44 -11.54 -36.45 -10.01 12.88
t4 0.00 0.00 17.11 41.90 60.58 32.76 68.31 58.97
t1 0.00 -0.00 -0.00 -0.00 -0.00 -32.95 -0.00 -32.72
Rente t2 -26.02 -28.40 -20.33 -24.27 -25.45 40.81 -27.05 -16.52
= 2% t3 14.74 37.36 -5.88 -10.13 -11.54 -36.45 -10.01 12.88
t4 0.00 0.00 17.11 41.90 60.58 32.76 68.31 58.97
Tabel 9.15: Payo matricer ved et volumen niveau på 0.0 MWh Ændringen i den optimale beslutningsregel påvirker kun simulationssti re. Ved en rente på 1% pumpes der gas ind til tidspunkt t1 , hvilket ikke sker ved en rente på 2%. Priserne tilbagediskonteres hårdere, hvilket ændrer estimaterne for β i den lineære regression for t1 . Disse er vist i tabel 9.16.
β0 β1 β2
Rente=1%
Rente=2%
-27.38 18.67 -0.83
-32.92 19.05 -0.83
Tabel 9.16 Ændringerne i β -estimaterne er ikke specielt store, hvilket heller ikke er forventet, da renten kun stiger med et procentpoint. Det huskes fra gennemgangen af LSMC, at ved et volumen niveau på 0% står valget mellem at pumpe den maksimale mængde gas ind i lageret, givet det ikke bryder den øvre rand, eller at gøre ingenting. Til at starte med udregnes continuation
86
value for renten på 1 % ved at forblive på volumen niveau på 0 %.
Cˆ 4 (t1 , 19.78607, 0) = −27, 38491 + 18.67161 · 19.78607 − 0.82870 · 19.786072 = 17.62510 (9.34) Hvis der pumpes ind i lageret ender volumen niveauet på 1.65 MWh. Dette er midt mellem de to gridpunkter 1.1 MWh og 2.2 MWh. Continuation value af at pumpe gas ind i lageret bliver derfor et vægtet gennemsnit af contination value for nabopunkterne. Disse er ikke udregnet her, men følger af fremgangsmåden fra ligning (9.34) og det illustrative eksempel.
Cˆ 4 (t1 , 19.78607, 0 + 1.65) = 40.39989 · 0.5 + 60.26957 · (1 − 0.5) = 50.33473
(9.35)
På baggrund af ovenstående continuation values er det nu muligt at træe det optimale valg.
Uinj (t1 , k = 1) = −19.78607 · 1.65 + 50.33473 = 17.68771
(9.36)
Unull (t1 , k = 1) = 17.62510
(9.37)
Ved en rente på 1% opnås det højeste payo ved at pumpe gas ind i lageret. De samme beregninger er gennemført for en rente på 2 %. Resultaterne er vist i tabel 9.17.
Cˆ 4 (t1 , 19.78607, 0 + 1.65) Cˆ 4 (t1 , 19.78607, 0)/Unull (t1 , k = 1) Uinj (t1 , k = 1)
Rente=1% 50.33473 17.62510 17.68771
Rente=2% 49.94506 17.36060 17.29804
Tabel 9.17 Det ses, at for en rente på 2%, kan det bedst betale sig ikke at pumpe gas ind i lageret. Desuden er de forventede payo generelt lavere ved en rente på 2 %, end ved en rente på 1 %. Dette skyldes den isolerede eekt af at skulle tilbagediskontere hårdere. Det er valgt at vise rentens ind ydelse på LSMC gennem denne simple gennemgang. Dette skyldes, at eekten af en renteændring kan være svær at analysere i et større setup, hvis der kun observeres priser på lageret. Skyldes prisændringen eekten gennem diskonteringen 87
isoleret, eller skyldes det også en ændring af den optimale beslutningsregel? Vi har i denne opgave valgt at benytte en rente på 0 %. Dette er et resultat af, at den faktiske rente, der kan observeres i markedet lige nu, er tæt på 0 %.
88
10
Sammenkobling af spotprismodellen og LSMC
Vi vil i dette afsnit beskrive sammenkoblingen af vores spotprismodel og LSMC metoden til prisfastsættelsen af gaslageret. Denne sammenkobling er ikke helt uproblematisk, da spotprismodellen afhænger af historikken i de forklarendevariable. I afsnit 5 fandt vi, at spotprisen afhænger af et 180 dages moving average af olieprisen, Degree Days, der fratrækkes gennemsnittet fra tidligere år, samt lagerniveauet, der ligeledes fratrækkes gennemsnittet fra de tidligere år. I LSMC metoden skal man i estimationen af continuation value på et givent tidspunkt t inkludere al tilgængelig information, der har relevans for den fremtidige prisudvikling. Dette medfører, at antallet af statevariable, som bør benyttes i estimationen, eskalerer. Betragt for eksempel situationen, hvor man ønsker at estimere continuation value for resten året den 1. april 2014. Spørgsmålet er,
hvilken tilgængelig information har ind y-
delse på prisudviklingen? For det første har 180 dages oliepriser, for det andet har de sidste 30 års temperaturer, som indgår i temperaturchokket og for det tredje har lagerniveauerne tilbage fra september 2007 en betydning. Oveni det kommer den autoregressive modellering af støjprocesserne. Grundet denne problematik vælger vi en pragmatisk tilgang, hvor der udvælges de tilstandsvariable, som vurderes at være de mest ind ydelsesrige. Vi har som udgangspunkt valgt at inkludere faktorerne moving average af olieprisen, temperaturchokket, lagerchokket og gasprisen på tid t. Denne grundmodel bliver speci ceret i begyndelsen af afsnit 11.
11
Resultater
Vi vil nu prisfastsætte et ktivt lager, teste forskellige aspekter og gennemgå hovedresultaterne. Lageret vi vil prisfastsætte har en maksimal volumen på 100.000 MWh, en injectionrate på 1.500 MWh, withdrawalrate på 3.000 MWh og renten antages at være 0. Afsnit 11.1 vil omhandle valget af basisfunktion, mens vi i afsnit 11.2 vil kigge på, hvor mange simulationer, der skal benyttes til prisfastsættelsen. Dernæst vil vi i afsnit 11.3 teste, hvordan opbygningen af volumegridet påvirker prisen. Dette sker ved at variere det brugervalgte hel89
tal precision-variablen. I afsnit 11.4 vil vi afprøve, hvordan LSMC metoden præsterer, når der benyttes forskellige tilstandsvariable. I afsnit 11.5 vil vi beskrive og undersøge størrelsen på et muligt bias i vores prissætning. I afsnit 11.6 vil vi kigge på, hvordan de realiserede prisstier fra 2012-2015 faktisk præsterer ved brug af LSMC metoden. Afslutningsvis vil vi sammenligne rene spotstrategier med handelsstrategier, der starter i et statisk hedge. Vi vil herunder angive Value at Risk og Expected Shortfall for de to strategier. I de forskellige analyser vil vi benytte en grundmodel. Denne består af gasprisen, moving average af olieprisen, temperaturchokket og lagerchokket til tid t, hvor samtlige kryds til og med anden orden er benyttet. Det giver 15 parametre. Desuden benyttes basisfunktionen monomials med 10.000 simulerede stier og precision på 1. Hvis ikke andet bliver nævnt i de følgende afsnit, så er grundmodellens input anvendt.
11.1 Basisfunktioner Dette afsnit omhandler valget af polynomium. Valget står mellem Monomials, Weighted Laguerre, Laguerre og Legendre. Det huskes fra afsnit 9.4.1, at der ved et større antal simulationer forventes ens resultater for Monomials, Laguerre og Legendre. For disse tre typer af polynomier er målet for nedenstående analyse ikke at vise grænsetilfældet, hvor polynomierne giver forskellige estimater. Det er derimod at vise, at ens resultater opnås, ved et stort antal af simulationer. Det store antal af simulationer sikrer, at der ikke er numeriske estimationsfejl. Af den grund benyttes 10.000 simulationer. I analysen benyttes gasprisen som eneste tilstandsvariabel. Det veri ceres, at ens resultater opnås for de forskellige typer af polynomier, hvor ordenen af polynomierne varieres. Ordenen af polynomierne går fra et 1. ordens til et 7. ordens polynomium. I tabel 11.1 ses resultaterne af simulationerne in-sample. In-sample angiver værdiansættelsen på baggrund af simulationerne. Udtrykket out-of-sample vil blive benyttet om den værdi, der opnås på den reelle realiserede prissti, når beslutningsreglen genereret af LSMC benyttes.
90
Orden Monomials W. Laguerre Laguerre Legendre
1
2
3
4
5
6
7
430.171 56.444 430.171 430.171
433.505 65.036 433.505 433.505
443.169 85.909 443.169 443.169
451.575 112.317 451.575 451.575
448.681 147.630 448.681 448.681
448.749 199.211 448.749 448.749
453.429 188.809 453.429 453.429
Tabel 11.1: Værdiansættelsen in-sample for lageråret 2014-2015 for varierende orden og type af polynomium. Monomials, Laguerre og Legendre polynomierne giver ens resultater, mens Weighted Laguerre in-sample giver et lavere payo på lageret. Weighted Laguerre performer relativt dårligt, hvilket kan overraske, da Boogert and de Jong [2011] får resultater, hvor denne performer næsten lige så godt som Monomials og Legendre. Det skal dog nævnes, at Boogert and de Jong [2011] benytter en anden prismodel, og at simulationerne og resultaterne er modelafhængige. Weighted Laguerre er i dette tilfælde ikke fordelagtig at vælge. Blandt de andre typer af polynomium kan der vælges frit. I opgaven vil vi benytte Monomials.
11.2 Antallet af simulerede pristier I dette afsnit vil vi teste, hvilken eekt antallet af simulationer har på prisfastsættelsen generelt. Hypotesen er, at prissætningen stabiliseres, når antallet af simulerede prisstier øges, da påvirkningen fra ekstreme priser mindskes. I tabel 11.2 stiger antallet af simulationer fra 2.500 til 28.500. I tabellen ses prisen og ændringen i prisen, når antallet af simulationer øges. Antal simulationer Pris Ændring
2.500 1.050.704
5.000 1.064.635 13.930
7.500 1.052.847 -11.787
10.000 1.053.886 1.039
12.500 1.061.547 7.661
15.000 1.060.588 -959
Antal simulationer Pris Ændring
17.500 1.061.078 490
20.000 1.056.981 -4.097
22.500 1.058.381 1.400
25.000 1.060.893 2.513
27.500 1.058.620 -2.273
28.500 1.059.419 799
Tabel 11.2: In-sample prisen givet antallet af simulationer. I tabellen ses det, at ændringen numerisk bliver mindre, når antallet af simulationer øges. Prisen svinger omkring 1.059.000, hvilket kan tyde på, at prisen konvergerer. Beregningstiderne for resultaterne i tabel 11.2 stiger lineært med antallet af simulationer. Prissætningen 91
ved 10.000 simulationer tog 83 minutter, mens det ved 28.500 simulationer tog 238 minutter. Denne tidsfaktor tages der højde for i de kontroller og tests, der foretages på LSMC algoritmen. Her vil benyttes 10.000 simulationer for at mindske beregningstiderne.
11.3 Precision I dette afsnit, vil det undersøges, hvordan opbygningen af volumegridet påvirker prisen. I LSMC algoritmen bestemmes afstanden mellem volumen gridpunkterne af forholdet mellem withdrawalraten og det brugervalgte heltal precision-variablen: wmax /precision. wmax er antaget konstant og eksogent bestemt. Det medfører, at volumegridet bliver mere nmasket, når precision-variablen øges. Det huskes fra afsnit 9.4, at når man ved en given volumeændring, opnår et volumeniveau uden for gridet, så bestemmes continuation value ved lineær interpolation af continuation values i nabopunkterne. Derfor bør et mere ntmasket net af punkter resultere i en skarpere prissætning af lageret. I prissætningen af lageret er injectionraten sat til 1.500 MWh/dag og withdrawalraten til 3.000 MWh/dag. Forholdet mellem de to rater har ligeledes betydning for prissætningen af lageret. Dette skyldes, at i langt størstedelen af gridpunkterne, vil det optimale valg stå mellem at pumpe den maksimaltilladelige mængde ind i lageret, gøre intet eller pumpe den maksimaltilladelige mængde ud af lageret. I tilfældet, hvor precision-variablen er sat til 4, er afstanden mellem punkterne i gridet wmax /precision = 3.000/4 = 750. Dette medfører, at man ved injection på 1.500 stort set altid vil ramme et punkt i gridet. Sænkes precision-variablen derimod til 3, vil der rammes uden for gridet ved max injection, da afstanden mellem volumen gridpunkterne i dette tilfælde er 1.000. I prissætningen af det ktive lager, vil det altså være fordelagtigt at vælge precision-variablen, så den går op i 2. Dog fremgår det af tabel 11.3, at precisionvariablen helst skal være over 2, da værdien af lageret her kun er 1.039.857. Kigges der på en precision over 2, så ses det, at prisen på lageret er højest, når precision er et lige heltal, og gridet rammes.
92
Precision Pris
1 1.053.886
2 1.039.857
3 1.060.380
4 1.062.728
5 1.061.154
Precision Pris
6 1.062.532
7 1.061.447
8 1.063.006
9 1.061.316
10 1.062.648
Tabel 11.3: In-sample prisen ved varierende precision. Det optimale valg af precision ændrer sig med setuppet for lageret. Flere faktorer skal overvejes. Udover forholdet mellem maksimal injection og withdrawal, tager vi ligeledes højde for beregningstider. En forhøjet precision-variabel øger stort set beregningstiden lineært. Ved en precision på 4 tog det 210 minutter at køre LSMC algoritmen. For en precision på 10 var beregningstiden steget til 496 minutter. I tabel 11.3 ses det, at prisen på lageret ved samme ændring i precision, falder med 80 EUR. Det er en meget lille ændring. Prisændringen sammenholdt med beregningstiden medfører, at precision på 4 virker som det optimale valg. I de resterende kontroller og test vil vi dog nøjes med at benytte en precision på 1 for at holde beregningstiden nede.
11.4 Tilstandsvariable I dette afsnit, vil vi undersøge, hvordan valget af tilstandsvariable i LSMC algoritmen påvirker prissætningen af lageret. Prissætningen er foretaget på de tre år gående fra 2012-2015. Her skal man være opmærksom på, at simulationerne benyttet til prissætning, er dannet ud fra prismodellen beskrevet i afsnit 5. Modellen baserer sig blandt andet på de realiserede stier fra 2012-2014, hvilket kan resultere i, at prissætningen i den periode har et positivt bias. I afsnit 7.1 er det vist, at modeller estimeret på få år har en høj RMSE. Den høje prædiktions fejl gør, at det er foretrukket, at prissætningen for årene 2012-2013 og 2013-2014 er fortaget in-sample. In-sample prissætningen giver muligvis værdiansættelsen et positivt bias, hvilket vil sige, at prisestimatet er lidt højere, end det vil være out-of-sample. Spørgsmålet er nu, hvilke tilstandsvariable skal benyttes for at optimere beslutningsreglen? Gasprisen er som udgangspunkt den mest oplagte tilstandsvariabel. Det testes derfor, hvordan gasprisen af varierende orden påvirker prisen på lageret. Resultaterne er vist i tabel 93
11.4.
PQ
q=1
Ptq
, Q :=
Pris in-sample 12-13 Pris in-sample 13-14 Pris in-sample 14-15
1
2
3
4
5
6
7
430.171 289.655 489.530
433.505 303.742 485.478
443.169 312.393 493.188
451.575 318.312 497.409
448.681 319.028 497.957
448.749 320.929 499.777
453.429 325.454 501.343
Tabel 11.4: Orden af monomials testet in-sample på gasårerne 2012-2015. Det ses, at prisen in-sample stiger, når ordenen øges. Dette er isoleret set et argument for at vælge en orden på syv eller højere. Det kan tænkes, at den støtte forbedring skyldes, at vi har benyttet de simulerede stier til både at danne beslutningsreglen, men også til udregningen af fremtidige payos. Den vedvarende forbedring kan derfor skyldes, at beslutningsreglen over tter til data. Dette problem diskuteres i afsnit 11.5. En anden problematik er, at for LSMC algoritmen kan et stort antal parametre give multikollinearitet i designmatricen. For at imødegå begge problematikker ønskes det at holde ordenen nede. Kigges på de tre år, bemærkes det, at den ekstra gevinst ved at øge ordenen, ser ud til at aftage. Ved en orden på 4 ser prisen ud til at have ramt et stabilt niveau. Næste skridt er at teste, hvordan tilstandsvariablene gaspris, oliepris, temperatur og lagerniveau i sammenspil påvirker lagerprisen. Som beskrevet i afsnit 10, afhænger vores spotprismodel af historikken i de forklarende variable. Derfor vil vi prisfastsætte lageret, hvor forskellige kombinationer tilstandsvariablene til tid t og t − 1 indgår i regressionen. Vi har valgt ikke at benytte ere laggede værdier af tilstandsvariablene for at holde antallet af parametre nede. For at forsimple notationen ses nedenfor en oversigt over de kombinationer af tilstandsvariable vi vil teste, hvor lag i N avnlag,o angiver om laggede værdier ingår i regressionen og o angiver ordenen. lag = 0 betyder derfor, at der ikke benyttes laggede værdier, mens lag = 1 angiver, at de laggede værdier fra tidspunktet t − 1 indgår ligeledes i regressionen.
• Alle0,1 : Beskriver en model, hvor tilstandsvariablene gaspris, oliepris, temperatur og lagerniveau indgår af første orden til tid t, hvilket giver en model med 5 parametre.
• U denGas0,1 , U denOlie0,1 , U denT emp0,1 , U denLager0,1 : Angiver re modeller, hvor fx 94
U denGas0,1 betyder, at olieprisen, temperaturen og lagerniveauet til tid t indgår som tilstandsvariable i regressionen af første orden. Hver af de re modellen indeholder 4 parametre.
• Alle0,2 : Angiver grundmodellen, hvor tilstandsvariablene indgår for tid t op til anden orden. Modellen indeholder 15 parametre.
• Alle1,1 : Denne model består af de re tilstandsvariable til tid t og t − 1. De indgår af første orden, hvilket giver 9 parametre.
• Alle1,2 : Angiver den fulde model med de re tilstandsvariable til tid t og t − 1 op til anden orden. Det giver 45 parametre.
• U denLager1,2 : Angiver modellen med gasprisen, olieprisen og temperaturen til tid t og t − 1 op til anden orden. Det giver 28 parametre. I tabel 11.5 og 11.6 ses resultaterne for 9 udvalgte modeller. Tilstandsvariable Pris 12-13 Pris 13-14 Pris 14-15
Alle0,1
U denT emp0,1
U denLager0,1
U denOlie0,1
U denGas0,1
1.093.763 890.860 1.032.625
1.068.244 796.579 1.011.800
1.067.584 888.312 965.073
669.757 539.053 704.795
1.004.498 711.289 948.668
Tabel 11.5: Antallet af tilstandsvariable testet in-sample på gasårerne 2012-2015.
Tilstandsvariable Pris 12-13 Pris 13-14 Pris 14-15
Alle0,2
Alle1,1
Alle1,2
U denLager1,2
1.107.058 909.519 1.053.886
1.372.341 1.073.991 1.192.531
1.395.291 1.096.356 1.204.938
1.360.671 1.090.450 1.108.073
Tabel 11.6: Antallet af tilstandsvariable testet in-sample på gasårerne 2012-2015. Det første, der bemærkes, er, at de 9 modeller for alle år prissætter lageret langt højere end modellerne kun indeholdende gasprisen af varierende orden fra tabel 11.4. Dernæst ses det, at værdien af lageret in-sampel generelt stiger, når der tilføjes ekstra parametre. Sammenlignes Alle1,1 og U denLager1,2 , ses det dog, at der er ere parametre i U denLager1,2 , men at 95
prisen er lavere for årene 12-13 og 14-15. Det tyder på, at det er informationen omkring lagerniveauet, der forbedrer in-sample prisfastsættelsen. Sammenlignes Alle0,1 og U denOlie0,1 fra tabel 11.5, ses det ligeledes, at olieprisen forbedrer in-sample prisen væsentligt. Vores grundmodel priser lageret til 1.053.886 EUR. Tilføjes viden omkring de laggede værdier af de re input, stiger værdien af lageret op til 1.204.938 EUR. Det er en væsentlig tilvækst i prisen, dog øges antallet af parametre fra 15 til 45. I de fortsatte tests vil vi stadig benytte grundmodellen, men til den endelig prisfastsættelse vil vi benytte Alle1,2 , da de laggede værdier forbedrer prissætningen.
11.5 Bias i prissætningen I LSMC algoritmen benyttes de simulerede prisstier til to ting. For det første estimeres
β -koe cienterne gennem en lineær regression på de simulerede prisstier. Dernæst træes beslutningerne på baggrund af continuation values. Hvis continuation value beregnes ud fra de samme stier, som er benyttet til at estimere β -koe cienterne, kan prisfastsættelsen have et positivt bias, hvilket vil sige, at den giver et højt prisestimat. Dette er nævnt i Stentoft [2004]. Spørgsmålene melder sig: Hvor stort er biaset? Hvordan påvirker valget af tilstandsvariable og antallet af simulationer biaset? Spørgsmålene adresseres ved at danne to grupper af simulerede prisstier. Den første gruppe vil benyttes til at estimere β -koe cienterne. Derefter anvendes beslutningsreglen på begge grupper separat. Biaset er så forskellen mellem den estimerede pris i hver af de to grupper. I afsnittet vil biaset analyseres gennem ændring af tilstandsvariablene og antallet af simulerede stier. I gur 11.1 er biaset plottet for 9 forskellige speci kationer. De to første speci kationer indeholder Alle0,2 og Alle0,1 , mens de resterende 7 speci kationer alene benytter gasprisen af varierende orden. Gruppe 1 og 2 består begge af 10.000 simulationer. Hypotesen er, at biaset bliver større, når antallet af parametre øges. Gruppe 1 vil altså over tte data og derfor estimere en højere pris end gruppe 2. Denne eekt ses tildels for gasprisen af første og anden orden, der har et lavere bias end de resterende modeller. Modellen Alle0,1 har det højeste bias på -7.754. Det huskes fra tidligere, at prisen på lageret 96
ligger omkring eller over 500.000 EUR. Set i det perspektiv er biaset altså relativt småt.
●
5000
●
●
●
P4t
P5t
●
●
0 −5000
EUR
●
●
●
Alle0.2
Alle0.1
P1t
P2t
P3t
P6t
P7t
Figur 11.1: Biaset ved 9 forskellige speci kationer, hvor Q i PtQ angiver den maksimale orden i polynomiet. Der er benyttet 2x10.000 simulationer. Antallet af simulationer påvirker ligeledes biaset. Hypotesen er her, at biaset nærmer sig nul, når antallet af simulationer øges. Det skyldes, at eekten af at bruge den enkelte prissti til begge ting, udvandes. Antallet af simulationer, som β -koe cienterne estimeres på baggrund af, øges fra 2.500 til 12.500. Samtidig indeholder gruppe 2 ved alle prissætninger de samme 10.000 observationer. I gur 11.2 ses biaset plottet, når antallet af simulationer øges.
97
8000 2000
●
●
−6000
−4000
−2000
0
Bias i EUR
4000
6000
●
●
2.500
5.000
7.500
●
10.000
12.500
Antal simulationer
Figur 11.2: Biaset ved et varierende antal simulationer. Tilstandsvariablene gaspris, oliepris, temperatur og lagerniveau af 1. orden er benyttet. Den numerisk laveste værdi af biaset ndes ved 12.500 simulationer, hvilket passer nt med hypotesen. På intervallet fra 5.000 til 12.500 simulationer ses ligeledes en tendens til, at biaset kommer tættere på nul. Biaset ved 2.500 simulationer falder lidt uden for hypotesen ved at have det numerisk andet laveste bias på -2.222. Alt i alt er der kun at sige, at biaset generelt ikke er specielt stort. Vi vil derfor fortsætte med at prisfastsætte simulationerne in-sample.
11.6 Resultater for de realiserede stier I de foregående afsnit blev der foretaget en modeltest in-sample. I dette afsnit vil der kigges på modellens præstation på de realiserede prisstier. Det vil sige, at der ses på, hvad der rent faktisk ville ske, hvis man følger en given strategi. Prissætningen er foretaget for årene gående fra 2012-2015. Det vil gøre det muligt at efterrationalisere på vores modelvalg. En god model bør vise de samme tendenser både in-sample og out-of-sample. Det skal dog nævnes, at out-of-sample værdien af lageret udelukkende baserer sig på en sti, nemlig den realiserede. Den realiserede prissti kan principielt helt bryde med de sammenhænge, som er 98
inkorporeret i modellen. Det kan den enten gøre, hvis der reelt set aldrig har været en eekt, eller hvis eekten ændrer sig. Denne problemstilling diskuterede vi i afsnit 7.1.1, speci kt eekten fra det store fald i olieprisen over gasåret 2014. Faldet var blandt andet et resultat af, at oliekartellet OPEC ville udkonkurrere den amerikanske skiferolie. I beslutningsreglen for den realiserede sti benyttes den realiserede oliepris som input. Dette er et eksempel på en eekt, som ændrer sig fra estimationsperioden til prisningsperioden. Forskellige valg af tilstandsvariable vil derfor præstere forskelligt med hensyn til både værdiansættelse in-sample og performance out-of-sample. I tabel 11.7 ses out-of-sample prisen, når den maksimale orden øges for et polynomium med gasprisen som tilstandsvariabel.
PQ
q=1
Ptq
1
2
3
4
5
6
7
356.933 -185.983 298.092
295.610 -282.068 222.268
295.936 -263.287 258.119
371.060 -320.183 328.535
335.820 -297.430 317.400
352.601 -315.468 295.126
311.003 -314.416 280.429
, Q :=
Out-of-sample 12-13 Out-of-sample 13-14 Out-of-sample 14-15
Tabel 11.7: Den maksimale orden af monomials polynomiet, som testes out-of-sample på gasårerne 2012-2015. I afsnit 11.4 blev det vurderet, at gasprisen af 4. orden virker som et fornuftigt valg, hvis der udelukkende vælges blandt polynomier med gasprisen som tilstandsvariabel. For årerne 12-13 og 14-15 havde denne model været det optimale valg, mens det i året 13-14 havde været det dårligste valg. Bedømt ud fra de tre out-of-sample stier, burde valget være faldet på modellen med gasprisen af 1. orden. I tabel 11.7 præsterer den bedst for året 13-14 og anden bedst for de sidste to år. Modellen med gasprisen af 4. orden blev in-sample overgået af modellerne med alle re tilstandsvariable. Out-of-sample prisen for disse er vist i tabellerne 11.8 og 11.9. Tilstandsvariable Pris out-of-sample 12-13 Pris out-of-sample 13-14 Pris out-of-sample 14-15
Alle0,1
U denT emp0,1
U denLager0,1
U denOlie0,1
U denGas0,1
482.966 -131.814 252.300
518.031 -272.627 237.354
424.727 -130.056 312.812
226.828 44.263 227.801
561.210 -272.859 273.943
Tabel 11.8: Out-of-sample prisen, når der benyttes forskellige tilstandsvariable.
99
Tilstandsvariable Pris 12-13 Pris 13-14 Pris 14-15
Alle0,2
Alle1,1
Alle1,2
U denLager1,2
455.362 -147.698 243.665
518.901 -128.934 252.091
536.062 -130.222 262.712
484.280 -144.078 331.019
Tabel 11.9: Antallet af tilstandsvariable testet out-of-sample på gasårerne 2012-2015. I tabel 11.8 ses det, at grundmodellen Alle0,2 præsterer middelmådigt i forhold til modellerne af 1. orden i tabel 11.8. For året 12-13 performer modellen uden gasprisen bedst med en pris på 561.210 EUR, mens for året 13-14 opnås den højeste pris ved at benytte modellen uden olieprisen. I 14-15 bliver prisen optimeret ved at fjerne lagerniveauet. Det må konstateres, at det har en betydning, hvilken model man vælger, men ligesom i in-sample testen, ser det ud til, at modellerne med laggede værdier tilfører værdi til lageret. For samtlige år præsterer modellerne Alle1,1 , Alle1,2 og U denLager1,2 bedre end vores grundmodel Alle0,2 . Dvs. inddragelsen af laggede værdier er vigtigt. Out-of-sample resultaterne ved at øge antal simulationer, valg af precision og valget af polynomium type er vedlagt i bilag under afsnittet 13.3.
100
11.7 Hedge og risikomål 11.7.1
Statisk hedge ved intrinsic-strategien
I indledningen til opgaven omtalte vi intrinsic-strategien, som baserer sig på at handle på sæsonudsvingene i forwardkurven. Alt gassen købes og sælges på en given dag og dette medfører, at strategien er riskofri, da man både kender købsprisen og salgsprisen. Derved er alle bevægelser ind og ud af lageret for det kommende år bestemt på handelsdagen. I denne opgave prises lageret over året gående fra den 1. april 2014 til den 31. marts 2015. Intrinsic-strategien bliver derfor vurderet ud fra forwardkurven den 31. marts 2014. Denne er plottet i gur 8.1. Strategien er meget intuitiv, da den i alt sin enkelthed går ud på at købe billigt og sælge dyrt. Rent praktisk kan strategien kræve en del kalkulationer grundet de tekniske begrænsninger for lageret. I Breslin et al. [2008] ndes den optimale intrinsicstrategi ved at opstille et maksimeringsproblem. Den optimale intrinsic-strategi ndes ved at løse nedenstående maksimeringsproblem.
XX
vlj Flj
(11.1)
vlj ≥ 0 X Il = vlj ≤ Imax
(11.2)
max vlj
j
l
Under bibetingelse af:
(11.3)
j
Wj =
X
vlj ≤ Wmax
(11.4)
l
Vl ≤ v max
(11.5) (11.6)
j>l
Hvor de nitionen af variablene er som følger:
• v max angiver den totale lagerkapacitet.
101
• Imax og Wmax er henholdsvis den maksimale månedlige injectionrate og withdrawalrate.
• Flj angiver spreadet mellem injection i måned l og withdrawal i måned j . • vlj er den købte mængde i spreadet. • Il angiver den mængde gas, der pumpes ind i lageret i måned l. • Wj angiver den mængde gas, der pumpes ud af lageret i måned j . • Vl er lagerniveauet i måned l. Bibetingelse (11.2) medfører, at der kun kan tages en lang position i et givent spread mellem to forwardkontrakter. Bibetingelserne (11.3) og (11.4) sikrer, at de injection og withdrawal positioner, der tages over en måned, ikke overstiger den maksimale månedlige injectionrate og withdrawalrate. Bibetingelsen (11.5) sørger for, at lagerkapaciteten i en given måned ikke overstiger den maksimale lagerkapacitet.
Ovenstående maksimeringsproblem vil nu anvendes på vores ktive lager, der har en størrelse på 100.000 MWh, en injectionrate på 1.500 MWh og withdrawalrate på 3.000 MWh. Det tager minimum 67 dage at fylde lageret op og 34 dage at tømme det igen. Der antages en rente på 0% og ingen transaktionsomkostninger. Grundet overlappende kontraktperioder benyttes intrinsic-strategien på forwardkontrakterne vist i tabel 11.10. I tabellen ses det, at lageret fyldes op over de tre første måneder. Lageret tømmes fuldt ud i Q1, da kvartalet har den højeste pris over året. Q1 indeholder 90 dage, hvilket er over de 34 dage, som er antallet af dage, der skal bruges til at tømme lagret.
102
Priser (EUR/MWh) Antal dage Volumen ændring (MWh) Volumen i lageret (MWh) Betalingsstrøm (EUR)
April
Maj
Juni
Q3
Q4
Q1
21,05 30 45.000 45.000 -947.250
21,25 31 46.500 91.500 -988.125
21,35 30 8.500 100.000 -181.475
21,60 92 0 100.000 0
24,02 92 0 100.000 0
24,93 90 -100.000 0 2.493.000
Tabel 11.10: Intrinsic strategien på lageråret 2014-2015. Nederst i tabel 11.10 ses betalingsstrømmene. Summeres disse fås et payo på 376.150 EUR. Dette payo kan sikres risikofri per den 31. marts 2014. Ulempen ved denne metode er, at bevægelserne i lageret i en vis udstrækning er låst. Derved har man mistet eksibiliteten til at handle, hvis man nder spotprisen fordelagtig. I ovenstående intrinsic-strategi bemærkes det, at der er 184 dage, hvor der ikke handles i lageret. Lageret er således fyldt ved udgangen af juni måned og tømmes i Q1, der dækker over månederne januar, februar og marts i 2015. Lageret skal være fyldt ved indgangen til januar, men spørgsmålet er, om det er muligt at få en gevinst ved at handle spotprisen i den mellemliggende periode? Det betyder, at gassen i lageret skal sælges og genkøbes i løbet af Q3 og Q4. Dette plejer ikke at være fordelagtigt, da gasprisen normalt er lavere i sensommeren end i vintermånederne, og man kæmper således mod sæsonudsvinget. For at teste, om der kan skabes pro t i disse måneder, er der foretaget modi kationer til LSMC koden. Denne sikrer, at der startes med et fyldt lager og sluttes med et fyldt lager. Den udglattede forwardkurve ses i gur 8.1, hvor det bemærkes, at prisen som forventet, stiger i perioden fra Q3 til Q4. I tabel 11.11 ses det forventede payo af lageret ved at handle på spotprisen i denne periode på 184 dage. Payo In-sample Out-of-sample
270.724 EUR 2.910 EUR
Tabel 11.11: Payo ved at handle i Q3 og Q4. In-sample resultatet er vurderet på 28.500 stier med speci kationen Alle1,2 . Payoet i tabel 11.11 er et gennemsnit af de enkelte stiers payo. Out-of-sample er LSCM-strategien kørt på den realiserede prissti for året gående fra 1. april 2014 til 31. marts 2015. Det ses, at der 103
er en markant forskel på in-sample og out-of-sample resultaterne. In-sample er vurderingen, at man kan tjene 270.724 EUR. Eekten af at handle i den mellemliggende periode viser sig dog ikke at øge værdien af lageret for den realiserede sti, da det kun giver en værdi på 2.910 EUR. Begge værdier skal ses i forhold til, at intrinsic-strategien sikrer et payo på 376.150 EUR oveni de henholdsvis 270.724 EUR in-sample og de 2.910 EUR out-of-sample.
11.7.2
Risikomål
Value at Risk Value at Risk (VaR) er et gængst benyttet risikomål. Value at Risk på et α-niveau angiver
1 − α fraktilet i pro t-loss fordelingen (P &L). Det indikerer således det tab, hvor sandsynligheden for at observere et større tab maksimalt må være 1 − α procent. I gur 11.3 ses et eksempel med en normalfordelt pro t-loss fordeling, hvor VaR på et 99% niveau er angivet. Dette svarer til 1% fraktilen i fordelingen. Value at Risk er givet ved ligning (11.7). VaRα (P &L) = sup{p&l | Prob(P &L ≤ p&l) ≤ 1 − α}
(11.7)
En kritik af VaR som risikomål er, at den blot angiver fraktilen, men ikke indeholder nogen information om fordelingen i den nedre hale i pro t-loss fordelingen. Således kan to pro tloss fordelinger med samme VaR have vidt forskellige fordelinger i halen. Et andet mål man kan benytte, og som indeholder nogen information om halens fordeling, er Expected Shortfall.
104
VaR
Figur 11.3: Value at Risk givet som det 1-α fraktil i pro t-loss fordelingen.
Expected Shortfall Expected Shortfall angiver det forventede tab, givet at tabet er blandt de 1 − α procent værste udfald. Dette kan skrives som i ligning (11.8).
ESα = E(P &L | P &L ≤ VaRα )
(11.8)
Denne betingede middelværdi vil i modsætning til VaR give en indikation om fordelingen i halen og dermed en fornemmelse af, hvor galt det kan gå, når det går galt. Med det sagt, så vil man stadig kunne have to forskellige udformninger af halen, der har det samme Expected Shortfall.
11.7.3
Risikomål benytte på LSMC værdiansættelsen
I den simulationsbaserede Least Squares Monte Carlo er 99% Value at Risk givet ved værdien af den sti, U b (t0 , Pt0 , v(t0 )), hvor der kun er 1% af de øvrige stier, som opnår et dårligere resultat i form af en lavere værdi. Expected shortfall på et 99% niveau er så gennemsnittet for de stier, der har en lavere værdi end den sti, som udgør 99% VaR. Vi vil udover disse to 105
størrelser også angive:
• Intrinsic værdien, som er værdien ved at følge intrinsic strategien. • Den realiserede værdi, som er den værdi, man ville opnå på den realiserede sti ved at følge beslutningsreglen, der er bestemt ud fra simulationerne.
• Værdien af den realiserede sti ved fuld information, hvilket er værdien, man vil kunne få ud af den realiserede sti, hvis man på starttidspunktet vidste, hvordan prisudviklingen vil blive resten af lageråret. Dette skal ses som en øvre grænse for lagerets værdi, der kan fås ved at handle spotprisen.
• Sandsynligheden for at opnå et tab er antallet af stier, der har en værdi under 0, divideret med det samlede antal stier.
• Sandsynligheden for at overstige intrinsic er antallet af stier, der har en højere værdi end intrinsic værdien, divideret med det totale antal simulerede stier.
• Sandsynligheden for at få en lavere værdi end den realiserede sti er antallet af stier, der har en lavere værdi end den realiserede værdi, divideret med det samlede antal stier.
11.7.4
Lagerårene 2012-2013 og 2013-2014
Der ses nu på værdiansættelserne og risikomålene for lagerårene 2012-2013 og 2013-2014. Dette gøres for at give et billede af, hvordan den realiserede værdi står mål med værdiansættelsen og de opstillede risikomål. For begge år er der brugt 10.000 simulationer. Opsætningen fra grundmodellen benyttes igen. Det vil sige, at tilstandsvariablene, der benyttes, er 180 dages moving average for olie, temperaturchokket og lagerchokket. De bruges op til 2. orden med monomials som basisfunktion. De 10.000 simulationer er foretaget ud fra modellen i ligning (5.6). Denne model er estimeret på baggrund af data i perioden fra den 17. september 2008 til 31. marts 2014, og værdiansættelsen skal derfor kategoriseres som en in-sample prisfastsættelse. Først ses på tilfældet, hvor der udelukkende handles på spotmarkedet. Dernæst 106
ses på tilfældet, hvor der på den sidste dag før lagerårets start laves et statisk hedge, som det blev beskrevet i afsnit 11.7.1. Der ndes med LSMC en handelsstrategi for de resterende dage, som ikke indgår i det statiske hedge. I tabel 11.12 ses, at for lageråret 2012-2013 er lagerets værdi for den rene spotstrategi estimeret til 1.107.058 EUR. Dette skal holdes op mod intrinsic værdien, der er på 433.500 EUR. Udregningen af intrinsic-strategien er vist i bilag i tabel 13.6. Værdien af lageret ved at følge den estimerede beslutningsregel på den faktiske prisudvikling er 455.361,7 EUR, hvilket næsten er 22.000 EUR mere end den risikofrie strategi. Til gengæld er den realiserede stis værdi markant lavere end værdiansættelsen på 1.107.058 EUR. 99% VaR angiver et tab på 401.886 EUR og ES et tab på 607.128. Dette indikerer, at man for de værste scenarier risikerer at inkassere markante tab. I tabel 11.12 er også opgivet, hvad den maksimale værdi af lageret er, hvis man kender hele den realiserede prissti på starttidspunktet. Dette svarer til, den værdi man kan få, hvis man er i besiddelse af en model, der er i stand til at forudsige stiens udvikling perfekt. Denne værdi er beregnet ved at løse det dynamiske problem beskrevet i afsnit 9.4, hvor man i stedet for at bruge den estimerede continuation value, bruger det tilbagediskonterede payo, som hver handling rent faktisk giver. Dette giver således en øvre grænse for, hvad den realiserede værdi kunne være blevet. Når der ses på den realiserede værdi og værdien af den realiserede sti med fuld information, skal det haves i mente, at det er bagudskuende information. Det kan benyttes til at evaluere sin model, men det er ikke information, som er kendt på dagen, hvor man fastsætter sin beslutningsregel for året. Det ses, at man på den realiserede sti maksimalt kunne have tjent 811.634 EUR, og at man derfor ikke ville kunne opnå værdien på de 1.107.058 EUR.
107
Nøgletal for lageråret 2012-2013 Værdiansættelse Intrinsic Value Realiseret værdi 99% VaR 99% ES Værdi af realiseret sti med fuld information Sandsynlighed for et tab Sandsynlighed for at overstige intrinsic Sandsynlighed for en lavere værdi end den realiserede
1.107.058,0 433.500,0 455.361,7 -401.886,8 -607.128,0 811.634,0 5,29 72,05 29,28
EUR EUR EUR EUR EUR EUR % % %
Tabel 11.12: Nøgletal for værdiansættelsen i lageråret 2012-2013. I tabel 11.13 ses resultaterne for 12-13, når der ved begyndelsen af året foretages et statisk hedge, som kun efterlader Q3 og Q4 åben for handel på spotprisen. Dette ses at have to konsekvenser. Den første er, at værdiansættelsen af lageret falder. Den anden er, at VaR og ES forbedres, hvilket indikerer, at risikoen mindskes. Dette ses også i gurerne 11.4 og 11.5. Værdien af den realiserede sti stiger dog med omtrent 14.000 EUR i forhold til den rene spot strategi. Sandsynligheden for et tab falder til 1,34 %. Denne størrelse var 5.29% uden hedget. Sandsynligheden for at overstige intrinsic værdien stiger, når man benytter hedget, fra 72,05 % til 92,11 %, men størrelsen på gevinsterne, som overstiger intrinsic bliver mindre jævnfør gur 11.5. Nøgletal for lageråret 2012-2013 med et statisk hedge Værdiansættelse Intrinsic Value Realiseret værdi 99% VaR 99% ES Sandsynlighed for et tab Sandsynlighed for at overstige intrinsic Sandsynlighed for en lavere værdi end den realiserede
737.774,0 433.500,0 477.142,5 -82.784,3 -428.044,3 1,34 92,11 10,46
EUR EUR EUR EUR EUR % % %
Tabel 11.13: Nøgletal for værdiansættelsen i lageråret 2012-2013 med udgangspunkt i et statisk hedge.
108
1.0
Expected Shortfall Value at Risk
0.0
0.2
0.4
Tæthed
0.6
0.8
Intrinsic value
−1.5 −1.0 −0.5
0.0
0.5
1.0
1.5
2.0
2.5
3.0
3.5
4.0
4.5
5.0
5.5
6.0
6.5
7.0
7.5
Millioner Euro
Figur 11.4: Fordelingen af lagerværdien for de 10.000 simulerede stier. Expected Shortfall og Value at Risk er udregnet på et 99% niveau.
Expected Shortfall 2.0
Value at Risk Intrinsic value
1.0 0.0
0.5
Tæthed
1.5
Uden hedge
−1.5 −1.0 −0.5
0.0
0.5
1.0
1.5
2.0
2.5
3.0
3.5
4.0
4.5
5.0
5.5
6.0
6.5
7.0
7.5
Millioner Euro
Figur 11.5: Fordelingen af lagerværdien for de 10.000 simulerede stier, når der startes med et statisk hedge. Expected Shortfall og Value at Risk er udregnet på et 99% niveau. 109
For lageråret 2013-2014 er billedet et noget andet end for 2012-2013. Her er spændet i forwardkurven meget snævert ved indgangen til lageråret. Dette skyldes dels, at gasprisen den 28. marts 2013 ligger på 34 EUR/MWh, hvilket må betragtes som et højt niveau for den del af året. Det gør, at det ikke er muligt at købe gassen billigt i starten af lageråret, hvilket mindsker sæsonspændet. I tabel 11.14 fremgår det, at intrinsic værdien er på 157.250 EUR, hvilket er en del lavere end for 2012-2013 året. VaR og ES viser tab på hhv. 831.242 EUR og 1.066.115 EUR. Nøgletal for lageråret 2013-2014 Værdiansættelse Intrinsic Value Realiseret værdi 99% VaR 99% ES Værdi af realiseret sti med fuld information Sandsynlighed for et tab Sandsynlighed for at overstige intrinsic Sandsynlighed for en lavere værdi end den realiserede
909.519,0 157.250,0 -147.698,1 -831.242,6 -1.066.115,0 367.985,0 15,03 75,71 10.09
EUR EUR EUR EUR EUR EUR % % %
Tabel 11.14: Nøgletal for værdiansættelsen i lageråret 2013-2014 . Værdiansættelsen bliver 909.519 EUR, hvilket er markant højere end intrinsic værdien. Udregningen af intrinsic-strategien er vist i bilag i tabel 13.7. På den realiserede sti inkasserer spotstrategien et tab på 147.698 EUR. Det er kun 10,09 % af de simulerede stier, som opnår et større tab. Ved fuld information omkring den realiserede værdi kan man maksimalt opnå en værdi på 367.985 EUR. Nøgletal for lageråret 2013-2014 med statisk hedge Værdiansættelse Intrinsic Value Realiseret værdi 99% VaR 99% ES Sandsynlighed for et tab Sandsynlighed for at overstige intrinsic Sandsynlighed for en lavere værdi end den realiserede
288.176,7 157.250,0 209.455,5 -467.936,7 -765.664,3 6,09 79,88 30,21
EUR EUR EUR EUR EUR % % %
Tabel 11.15: Nøgletal for værdiansættelsen i lageråret 2013-2014 med statisk hedge. 110
Tabel 11.15 indeholder nøgletal for lageråret 2013-2014, når der ved begyndelsen af året er foretaget et statisk hedge. Intrinsic-strategien, fra tabel 13.7, låser mulighederne for at handle i juni, Q3 og Q1. Der kan således kun handles på spotprisen i april, maj samt oktober, november og december. Hedget forbedrer også denne gang værdien af den realiserede sti, der går fra -147.698,1 EUR til 209.455,5 EUR. VaR og ES forbedres ligeledes. På baggrund af backtestene i årene 2012-2013 og 2013-2014 vil vi anbefale, at spotstrategien startes i et delvist eller fuldt hedge. Dette skyldes primært tabet på 147.698 EUR, som man ville have
1.0
inkasseret ved den rene spotstrategi, hvor der ikke startes i et hedge.
Expected Shortfall Value at Risk
0.0
0.2
0.4
Tæthed
0.6
0.8
Intrinsic value
−1.5
−0.5
0.0
0.5
1.0
1.5
2.0
2.5
3.0
3.5
4.0
4.5
5.0
5.5
6.0
6.5
7.0
7.5
Millioner Euro
Figur 11.6: Fordelingen af lagerværdien for de 10.000 simulerede stier. Expected Shortfall og Value at Risk er udregnet på et 99% niveau.
111
3.0
Expected Shortfall Value at Risk
2.5
Intrinsic value
1.5 0.0
0.5
1.0
Tæthed
2.0
Uden hedge
−1.5 −1.0 −0.5
0.0
0.5
1.0
1.5
2.0
2.5
3.0
3.5
4.0
4.5
5.0
5.5
6.0
6.5
7.0
7.5
Millioner Euro
Figur 11.7: Fordelingen af lagerværdien for de 10.000 simulerede stier, når der startes med et statisk hedge. Expected Shortfall og Value at Risk er udregnet på et 99% niveau.
11.7.5
Lageråret 2014-2015
Nu ses der på værdiansættelsen i lageråret 2014-2015. Vi vil på baggrund af analysen i afsnit 11.4 benytte opsætningen Alle1,2 , hvilket er alle tilstandsvariablene til og med et lag, der inddrages op til 2. orden. Begge strategier vurderes ud fra de samme 28.500 simulationer, og lageret har en injectionrate på 1.500 MWh, en widrawalrate på 3.000 MWh og en maksimal kapacitet på 100.000 MWh. I tabel 11.16 ses værdiansættelsen at være på 1.213.935 EUR, hvilket er mere end intrinsic værdien på 376.150 EUR. Man løber dog en betydelig risiko, hvilket ses af 99% Value of Risk og 99% Expected Shortfall estimaterne. VaR angiver, at man med en sandsynlighed på 1% vil realisere et tab, der er større end 19.068 EUR. Givet, at et af de 1% værste scenarier realiseres, forventes et tab på 109.942 EUR.
112
Nøgletal for lageråret 2014-2015 Værdiansættelse Intrinsic Value Realiseret værdi 99% VaR 99% ES Værdi af realiseret sti med fuld information Sandsynlighed for et tab Sandsynlighed for at overstige intrinsic Sandsynlighed for en lavere værdi end den realiserede
1.213.935,0 376.150,0 268.401,5 -19.068,4 -109.941,7 983.012,0 1,25 82,49 10,98
EUR EUR EUR EUR EUR EUR % % %
Tabel 11.16: Nøgletal for værdiansættelsen i lageråret 2014-2015. I gur 11.8 er fordelingen af de enkelte stiers værdi illustreret ved hjælp af et histogram og en udglattet tæthed, hvilket giver indtrykket af, at der er en betydelig variation i værdien ved at benytte en ren spotpris strategi. I tabel 11.16 er også angivet sandsynligheden for at opnå en negativ værdi, der er på 1,25%, samt sandsynligheden for at opnå en højere værdi end intrinsic-strategien, der på 82,49%. Det at benytte en ren spotpris strategi er ligesom i lagerårene 2012-2013 og 2013-2014 risikabelt. Vi vil nu sammenholde dette med tilfældet, hvor der startes med det statiske hedge, som vi har beskrevet i afsnit 11.7.1. For det statiske hedge er der 184 dage fra den 1. juli 2014 til den 31. december 2014, hvor man kan handle spotprisen. Betingelsen er her, at man starter med et volumeniveau på 100% den 1. juli, hvilket man også skal slutte med den 31. december. Dette skyldes, at man har forpligtet sig på forwardmarkedet ved at sælge hele lagerets kapacitet i perioden fra d. 1. januar 2015 til den 31. marts 2015. I tabel 11.17 er værdiansættelsen angivet, når der startes med et statisk hedge. Som i det tidligere tilfælde er 99% Value at Risk og Expected Shortfall angivet. Det vurderes, at man kun i 1% af tilfældene vil opnå et resultat, der er mindre end 85.044 EUR, og at det gennemsnitlige resultat for de 1% værste scenarier er et tab på -117.128 EUR. Så hvor VaR forbedres, forværres ES med 8.000 EUR. Dette kan skyldes, at der handles mod sæsonen i de 184 dage, som hedget ikke har låst, hvilket gør, at der med hedget er et worst case scenario på -2 millioner EUR, som skiller sig ud fra de andre. Dette er et større tab end worst case scenario uden hedge. Dog ses det i gur 11.9, at variansen falder i fordelingen af priser, når det statiske hedge indgås. 113
Nøgletal for lageråret 2014-2015 med et statisk hedge Værdiansættelse Intrinsic Value Realiseret værdi 99% VaR 99% ES Sandsynlighed for et tab Sandsynlighed for at overstige intrinsic Sandsynlighed for en lavere værdi end den realiserede
646.873,9 376.150,0 379.060,0 85.043,7 -117.127,9 0,62 93,24 6,95
EUR EUR EUR EUR EUR % % %
Tabel 11.17: Nøgletal for værdiansættelsen i lageråret 2014-2015, der startes i et statisk hedge. Sandsynligheden for at opnå et tab falder en smule fra 1,25 % til 0,62%. Sandsynligheden for at overgå intrinsic-strategien er nu forøget til 93,24 %. Ved den rene spotmarked strategi er sandsynligheden for at overstige intrinsic værdien altså lavere, men gevinsten er markant højere, når den overstiger intrinsic værdien. Dette kan ses ved at sammenholde gur 11.8 og 11.9.
0.6
Expected Shortfall Value at Risk
0.3 0.0
0.1
0.2
Tæthed
0.4
0.5
Intrinsic value
−0.5
0.0
0.5
1.0
1.5
2.0
2.5
3.0
3.5
4.0
4.5
5.0
5.5
6.0
6.5
Millioner Euro
Figur 11.8: Fordelingen af lagerværdien for de 28.500 simulerede stier. Expected Shortfall og Value at Risk er udregnet på et 99% niveau.
114
2.5
Expected Shortfall Value at Risk Intrinsic value
0.0
0.5
1.0
Tæthed
1.5
2.0
Uden hedge
−0.5
0.0
0.5
1.0
1.5
2.0
2.5
3.0
3.5
4.0
4.5
5.0
5.5
6.0
6.5
Millioner Euro
Figur 11.9: Fordelingen af lagerværdien for de 28.500 simulerede stier, når der startes med et statisk hedge. Expected Shortfall og Value at Risk er udregnet på et 99% niveau. De to tilfælde ligger i hver ende af skalaen med hensyn til risikoaversion, hvor den rene spotstrategi er meget risikovillig, er strategien startende i det statiske hedge risikoavers. En mellemting kan opnås, hvis man kun laver et statisk hedge på en delmængde af lagerets fulde kapacitet og så anvender en spotstrategi på den overskydende kapacitet. Dette vil dog kræve modi kationer i den eksisterende kode, da man ved at handle på forwardmarkedet på forhånd optager en del af de tilladte withdrawal og injection mængder. Man skal således med udgangspunkt i det statiske hedge nde ud af, hvor meget overskudskapacitet, der er tilbage for hver given dag. Dette gør, at de maksimale withdrawal- og injectionrater ikke er konstante over tid, hvilket programmet lige nu er baseret på.
11.7.6
Betydningen af reduktionen af stier
I afsnittet 7.3 blev det beskrevet, hvordan nogle af de simulerede stier bliver ekstreme, hvilket primært tilskrives enhedsroden i oliemodellen. Som alternativ til at opstille en ny model, motiverede vi at fjerne de henholdsvis 2.5 % største prisstigninger og prisfald hen 115
over året. Det er disse stier, der vil drive de største tab og gevinster. Beslutningsreglen vil forsøge at undgå de store tab og eksponere sig selv mod de store gevinster, i stedet for at tilpasse sig de mere realistiske gennemsnitlige stier, hvor olieprisen ikke styrtdykker eller er på himmel ugt. Det vil her blive illustreret, hvilken eekt det vil have, ikke at fjerne de ekstreme stier. Der ses igen på lageråret 2014-2015, og der benyttes det samme setup som i den tidligere værdiansættelse af dette lagerår, dog benyttes der kun 10.000 simulationer her. I guren 11.10 ses histogrammet og den udglattede tæthed, når der indgås det statiske hedge. I guren er også den udglattede tæthed for den ren spotstrategi, der nu har tungere
2.5
haler end tidligere.
Expected Shortfall Value at Risk Intrinsic value
0.0
0.5
1.0
Tæthed
1.5
2.0
Uden hedge
−2
0
2
4
6
8
Millioner Euro
Figur 11.10: Fordelingen af lagerværdien for de 10.000 simulerede stier, når der ikke fjernes ekstreme simulationer. Ved at sammenholde nøgletallene for værdiansættelserne med og uden de ekstreme situationer, ses det, at værdien generelt stiger, når man inkluderer de ekstreme situationer. Dette ses ved at sammenholde tabellerne 11.16 og 11.17 med 11.18 og 11.19. VaR og ES indikerer større tab, hvilket giver mening, da de laveste prisscenarier her er inkluderet. Den realiserede værdi out-of-sample stiger med 40.000 EUR i forhold til, når der fjernes stier. 116
Nøgletal for lageråret 2014-2015 med ekstreme simulationer Værdiansættelse Intrinsic Value Realiseret værdi 99% VaR 99% ES Sandsynlighed for et tab Sandsynlighed for at overstige intrinsic Sandsynlighed for en lavere værdi end den realiserede
1.270.527,0 376.150,0 307.281,4 -69.137,3 -156.880,7 2,13 77,51 18,00
EUR EUR EUR EUR EUR % % %
Tabel 11.18: Nøgletal for værdiansættelsen i lageråret 2014-2015 med ekstreme simulationer.
Nøgletal for lageråret 2014-2015 med ekstreme simulationer og hedge Værdiansættelse Intrinsic Value Realiseret værdi 99% VaR 99% ES Sandsynlighed for et tab Sandsynlighed for at overstige intrinsic Sandsynlighed for en lavere værdi end den realiserede
663.072,6 376.150,0 386.950,0 48.516,4 -184.396,4 0,79 92,51 8,26
EUR EUR EUR EUR EUR % % %
Tabel 11.19: Nøgletal for værdiansættelsen i lageråret 2014-2015 med ekstreme simulationer, der startes i det statiske hedge.
117
11.8 Opsamling af hovedresultater Vi fandt i afsnit 11.1, at basisfunktionerne Laguerre, Legendre og Monomials præsterer bedre end Weighted Laguerre. Blandt de tre basisfunktioner er vi indierente, da de prisfastsætter lageret ens. Vi viste i afsnit 11.2, at prisfastsættelsen afhænger af antallet af simulationer. Estimationen af lagerprisen bliver mere stabil, når der benyttes et større antal simulationer. Volumen gridet i LSMC metoden har ligeledes betydning for prissætningen. I afsnit 11.3 viste vi, at mindskes afstanden i volumen gridet, så øges prisen af lageret. Vi viste ligeledes, at prisfastsættelsen afhænger af forholdet mellem withdrawalraten og injectionraten, hvilket betyder, at afstanden mellem punkterne i volumen gridet skal tilpasses det speci kke lager, der prisfastsættes. I afsnit 11.4 testede vi LSMC metoden med forskellige tilstandsvariable. Analysen viste, at gasprisen, moving average af olieprisen, lagerchokket og temperaturchokket skal indgå i LSMC prissætningen. Historikken i de re tidsrækker forbedrede ligeledes prissætningen, hvilket vi viste, da vi indragede laggede værdier. I afsnit 11.5 kiggede vi på et muligt bias i prisfastsættelsen. Dette blev testet ved at variere tilstandsvariablene og antallet af simulationer. Biaset viste sig at have en meget lille betydning for prissætningen. I afsnittet 11.6 viste vi, hvordan de realiserede stier præsterer på åreren 2012-2015. Resultaterne viste her nogle af de samme tendenser som in-sample. Historikken i tilstandsvariablene havde for alle årerne betydning for prissætningen. Afslutningsvis fortog vi en analyse af risikopro len for LSMC strategien sammenlignet med intrinsic-strategien. Denne viste, at handel på spotprisen alene er meget risikofyldt. En blandet strategi virker som det fornuftige valg. En del af lageret reserveres til handel på forwardmarkedet og den resterende del handles på spotmarkedet. Derudover bemærkes det, at in-sample priserne generelt er ret høje. Værdien af lageret med fuld information er for alle årene lavere end den værdiansættelse, som fremkommer, ved at benytte LSMC metoden. Det er gældende for alle årene og er måske et tegn på, at den realiserede sti giver lavere priser end den udglattede forwardkurve. I afsnit 11.7.6 viste vi, hvordan prissætningen ville have været, hvis vi ikke havde reduceret antallet af simulerede stier. Den realiserede sti for året 2014-2015 præsterer i dette tilfælde bedre med de ekstreme værdier. På den ene side giver dette resultat anledning til at beholde 118
samtlige stier. På den anden side viste vi i afsnit 7.3, at gennemsnittet af 2,5 % højeste observationer er seks gange større end prisen den 31. marts. 2014. Konklusionen er, at man ikke bør reducere i antallet af stier, men i stedet overveje en ny model for olieprisen, som i vores modellering har en enhedsrod. Dette er en bedre tilgang end at fjerne prisstier, hvilket man egentlig ikke bør gøre, men som vi har gjort grundet de ekstreme oliestier. I analysen af LSMC metoden har vi antaget, at der ikke er nogle transaktionsomkostninger ved at handle. Vi har ligeledes negligeret bid-ask spreadet og eventuelle omkostninger ved at vende strømmen af gas fra at pumpe ind til at pumpe ud og omvendt. Disse antagelser er lavet for at kunne fokusere vores studier på andre dele prissætningen, men det er bestemt et område af opgaven, der kan videreudvikles på.
119
12
Konklusion
I denne afhandling understøtter vi empirisk, at gaslagerbeholdningen har en betydning for udviklingen i spotprisen. Vi fremsætter en spotprismodel, hvor vi med inspiration fra modellerne beskrevet i Stoll and Wiebauer [2010] og Müller et al. [2015], inkluderer gaslagerbeholdningen. Som i Müller et al. [2015] nder vi en positiv korrelation mellem olieprisen og gasprisen via et 180 dages moving average. Müller et al. [2015] tilskriver dels eekten, at gasprisen indekseres på olieprisen og dels, at olieprisen følger konjunkturerne, som også driver efterspørgslen efter gas. Vores out-of-sample test af spotprismodellen i lageråret 20142015 viser, at olieprisens styrtdyk ikke nødvendigvis leder til et tilsvarende dyk i gasprisen, når konjunkturerne er gunstige. Eekten af temperatur viser som forventet, at gasprisen er negativt korreleret med Degree Days. Lagerbeholdningen påvirker gasprisen således, at lagerniveauer, der ligger under gennemsnittet for de tidligere år, vil medføre en stigning i gasprisen. Eekten viser sig at være sæsonafhængig. Eekten er lavest i maj, hvor et fald i lagerbeholdningen med 1 procentpoint vil medføre en stigning på 0.8 procentpoint i gasprisen. Den største eekt ligger i overgangen mellem november og december, hvor et tilsvarende fald på 1 procentpoint vil give en stigning i gasprisen på lidt over 1.8 procentpoint. Dernæst opstiller vi tidsrækkemodeller for olieprisen, temperaturen og gaslagerbeholdningen, som indgår i simulationen af gaspriser. Disse simulationer kalibrerer vi til en udglattet forwardkurve, hvor vi i udglatningen følger Benth et al. [2007]. Kalibreringen bliver foretaget på den sidste dag inden lagerårets start. De kalibrerede prissimulationer anvendes til prisfastsættelsen af en gaslagerkontrakt ved Least Squares Monte Carlo metoden, som er beskrevet i Boogert and De Jong [2008]. Vi viser i den forbindelse, at man løber en betydelig risiko, hvis man ikke gør brug af et statisk hedge. Dette baserer sig på vurderingen af Value at Risk, Expected Shortfall samt fordelingen af lagerværdien. Handelsstrategien med det statiske hedge klarer sig generelt bedre end den rene spotmarkedsstrategi både in- and outof-sample. Vi har vist at den blandede strategi er pro tabel, da den opnår et bedre resultat end den rene intrinsic strategi både in- and out-of-sample.
120
13
Bilag
Dette afsnit vil indeholde udvalgte tabeller, gure og kodestykker. Vi har programmeret i
R, som kan downloades gratis via http://www.r-project.org/. En lang række af modi kationerne til LSMC koden er udeladt, men kan tilsendte på anmodning.
13.1 Prismodel
Intercept s(t) sin cos t N HDDt Ψt
Model Müller Estimate Std. Error -3.52 0.28 0.27 0.11 1.29 0.08 -0.19 0.08 -0.00 0.00 0.01 0.00 0.36 0.01
Model Stoll et al. Estimate Std. Error 10.87 0.27 0.26 0.20 1.33 0.13 0.62 0.13 0.01 0.00 0.01 0.00 -
Tabel 13.1: Estimater for modellerne af Müller og Stoll et al. Estimaterne til den skæve t-fordeling fra simulationsmodellen af olieprisen.
µ θ ω α1 β1 skew shape
Estimate 0.00 -0.05 0.00 0.04 0.96 0.94 9.33
Std. Error 0.00 0.02 0.00 0.01 0.01 0.02 1.29
t value 2.03 -3.02 2.07 7.26 171.26 44.20 7.23
Pr(>|t|) 0.04 0.00 0.04 0.00 0.00 0.00 0.00
Tabel 13.2: Parameter estimater for modellen i ligningssystemet (6.3).
121
22 20 18 5 0 70 65 0.00 0.05 0.10 −0.10
Lagerchok
60
EUR/bbl
75
80 −5
Degree Days
10
16
EUR/MWh
24
13.2 2014-2015
apr 14
maj 14
jun 14
jul 14
aug 14
sep 14
okt 14
nov 14
dec 14
jan 15
feb 15
mar 15
Figur 13.1: Realiserede værdier af gasprisen, temperaturchokket, 180 dages moving average af olieprisen og lagerchokket.
13.3 Resultater I tabel 13.3 ses, hvordan de forskellige basisfunktioner performer på den realiserede prissti for året 2014-2015. Orden Monomials W. Laguerre Laguerre Legendre
1
2
3
4
5
6
7
298.092 43.865 298.092 298.092
222.268 74.763 222.268 222.268
258.119 38.038 258.119 258.119
328.535 139.916 328.535 328.535
317.400 234.198 317.400 317.400
295.126 -76.083 295.126 295.126
280.429 -219.068 280.429 280.429
Tabel 13.3: Out-of-sample prisen ved varierende orden og type af polynomium. I tabel 13.4 ses prissætningen af lageret out of sample på året 2014-2015, når der ændres på precision-variablen.
122
Precision Pris
1 243.665
2 243.705
3 248.035
4 237.922
5 246.436
Precision Pris
6 232.850
7 246.015
8 243.114
9 256.453
10 249.956
Tabel 13.4: Out-of-sample prisen ved varierende precision. I tabel 13.5 ses out-of-sample prisen ved et varierende antal simulationer. Antal simulationer Pris
2.500 230.986
5.000 232.111
7.500 236.759
10.000 243.665
12.500 240.498
15.000 242.595
Antal simulationer Pris
17.500 245.409
20.000 243.599
22.500 246.009
25.000 247.326
27.500 249.827
28.500 249.144
Tabel 13.5: Out-of-sample prisen givet antallet af simulationer. Intrinsic-strategien for 2012-2013 ses i tabel 13.6.
Priser (EUR/MWh) Antal dage Volumen ændring (MWh) Volumen i lageret (MWh) Betalingsstrøm (EUR)
April
Maj
Juni
Q3
Q4
Q1
25,15 30 45.000 45.000 -1.131.750
25,45 31 46.500 91.500 -1.183.425
25,45 30 8.500 100.000 -216.325
25,7 92 0 100.000 0
28,6 92 0 100.000 0
29,65 90 -100.000 0 2.965.000
Tabel 13.6: Intrinsic-strategien for 2012-2013 Intrinsic-strategien for 2013-2014 ses i tabel 13.7.
Priser (EUR/MWh) Antal dage Volumen ændring (MWh) Volumen i lageret (MWh) Betalingsstrøm (EUR)
April
Maj
Juni
Q3
Q4
Q1
29,35 30 0 0 0
27,35 31 0 0 0
27,00 30 45.000 45.000 -1.215.000
27,05 92 55.000 100.000 -1.487.750
28,10 92 -0 100.000 0
28,60 90 -100.000 0 2.860.000
Tabel 13.7: Intrinsic-strategien for 2013-2014
13.4 Kode til prismodellen Spotpris modellering 123
1
library (" forecast "); library (" tseries ") ; library ( " R.utils ") ; library (" xtable ")
2 3
dat ← read.csv ( "F: / April 2014 Data / Samlet data APRIL 2014 .csv " , header = TRUE , sep =" ;" , dec = " ," )
4 5
dat1 ← dat [, c (3:7 ,16 ,17 ,21) ] ; NROW ( dat1 ) ; head ( dat1 )
6
dat2 ← dat [, c (3:7 ,16:19) ] ; head ( dat2 )
7
dat3 ← dat [, c (4:7 ,9 ,16 ,17 ,18 ,21) ] ; head ( dat3 )
8
dat4 ← dat [, c (3:7 ,16 ,18 ,21) ] ; head ( dat4 )
9 10
model2s ← lm ( Outlier.korrigeret.data∼ .-X.180.dages.MA..med.interpolation-JM_MED_CDD-W_shocks , data = dat2 , na.action = na.omit )
11
summary ( model2s )
12 13 14
# ### Model fit
15
model1 ← lm ( log ( Outlier.korrigeret.data )∼. +( Storage_andel_ugentlig ) *( Cos + Sin )+ log ( X.180.dages.MA..med.interpolation ) -X.180.dages.MA..med.interpolation , data = dat1 , na.action = na.omit )
16
summary ( model1 )
17
# saveObject ( model1 ,"/ Volumes / USB / Parameterestimater / Model1.txt ")
18
y ← resid ( model1 )
19 20
# #############
21
IC← matrix ( numeric () , nrow =7 , ncol =2)
22
rownames ( IC ) ← paste (" AR " ,1:7 , sep = "")
23
colnames ( IC ) ← c( " AIC " ," BIC ")
24
124
25
for ( i in 1:7) {
26
IC [i ,1]
← Arima (y ,c(i ,0 ,0) , include.mean = FALSE )$ aic
27
IC [i ,2]
← Arima (y ,c(i ,0 ,0) , include.mean = FALSE )$ bic
28
}
29
IC
30
xtable (IC , digits =0)
31
( AR5 ← Arima (y ,c (5 ,0 ,0) , include.mean = FALSE ))
32
res ← resid ( AR5 )
33
adf.test ( res )
34
# saveObject ( AR5 ,"/ Volumes / USB / Parameterestimater / Gas_AR5.txt " )
35 36
# ########### Tunge haler
37
# ## AIC mere stringent.
38
stepAIC ← stepAIC.ghyp ( res )
39
stepAIC $ fit.table
40 41
nig ← fit.NIGuv ( res )
42
x ← seq ( min ( res ) , max ( res ) , length.out =1000)
43 44
# VARIANCE GAMMA
45
hyp_vg ← ghyp ( lambda = vg@lambda , mu = vg@mu , sigma = vg@sigma , gamma = vg@gamma )
46
saveObject ( nig , "/ Volumes / USB / Parameterestimater / NIG.txt " )
47 48
modelStoll ← lm ( Outlier.korrigeret.data∼ .-X.180.dages.MA..med.interpolation-Storage_andel_ugentlig , data = dat4 )
49
summary ( modelStoll )
50
modelMuller ← lm ( Outlier.korrigeret.data∼.-Storage_andel_ugentlig , data = dat4 )
125
51
summary ( modelMuller )
Listing 1: Spotpris modellering
Temperaturmodel 1
dtemp ← read.csv ( file.choose () , header = TRUE , sep =" ;" , dec =" ,")
2 3
dtemp1 ← cbind ( dtemp , sin (2* pi * dtemp $ t_fra_1950 / 365 .25 ) , cos (2* pi * dtemp $ t_fra_1950 / 365 .25 ))
4
names ( dtemp1 ) ← c( " temp " ,"t " ," sin " ," cos " )
5
head ( dtemp1 )
6
fit ← lm ( temp∼. , data = dtemp1 )
7 8
( ARMA ← Arima ( dtemp1 [ ,1] , c (4 ,0 ,0) , xreg = dtemp1 [ , -1 ]) )
9
cbind ( ARMA $ coef , sqrt ( diag ( ARMA $ var.coef )) )
10 11
ARres ← resid ( ARMA )
12 13
# ################## Simulation af temperatur
14
ny_temp ← dtemp1 $t [ NROW ( dtemp1 )] + 1:365
15
ny_reg ← cbind ( ny_temp , sin (2* pi * ny_temp / 365 .25 ) , cos (2* pi * ny_temp / 365 .25 ))
16
colnames ( ny_reg ) ← c( "t" ," sin " ," cos ") ; head ( ny_reg )
17 18
# write.table ( ny_reg ," / Volumes / USB / Simulation / Ny_reg_temp.txt " )
19
AR4 ← list ( ARMA $ coef )
20
saveObject ( ARMA ,"/ Volumes / USB / Parameterestimater / Temp_AR4.txt " )
21 22
sim_temp ← simulate.Arima ( ARMA , nsim = NROW ( ny_reg ) , xreg = ny_reg )
23
dat14 ← read.csv ( file.choose () , header = TRUE , sep =" ;" , dec =" ," )
24
126
25
simuler_Temp ← function () {c ( ( dat14 [1 ,2] + cumsum ( pmax (15 -sim_temp ,0) + pmax ( sim_temp-21 ,0) ) - dat14 [ ,1]) [1:183] , ( cumsum ( pmax (15 -sim_temp [184:365] ,0) + pmax ( sim_temp [184:365] -21 ,0) ) - dat14 [184:365 ,1]) ) }
Listing 2: Temperaturmodel
Lagermodel 1
Storage_reg
← read.csv ( "/ Volumes / USB / April 2014 Data / Storage
sin-cos-reg.csv " , header = TRUE , sep =" ;" , dec =" ,") [1:2388 , -c (1:2 ,7:11) ] 2
dato ← read.csv ( "/ Volumes / USB / April 2014 Data / Storage sin-cos-reg.csv " , header = TRUE , sep =" ;" , dec =" ,") [1:2388 ,1]
3
pred_cos_sin ← tail ( read.csv ( "/ Volumes / USB / April 2014 Data / Storage sin-cos-reg.csv " , header = TRUE , sep =" ;" , dec =" ,")[ , -c (1:2 ,7:11) ] ,371)
4
pred_cos_sin ← pred_cos_sin [ seq (7 ,365 ,7) ,-5 ]
5 6
u.storage.reg ← Storage_reg [ seq (1 ,2388 ,7) ,]
7
attach ( u.storage.reg )
8 9 10
plot.ts ( Storage_andel ) St_reg ← u.storage.reg [, - (5:6) ] ; head ( St_reg ); tail ( St_reg )
11 12
# ################## Logit transformation
13
logit ← log ( Storage.procent.ugentlig / (1 -Storage.procent.ugentlig ))
14 15
# ################# Model NY
16
library ( R.utils )
17
ARMA00 ← Arima ( logit ,c (0 ,0 ,0) , seasonal = list ( order = c (0 ,0 ,0) , period =7) , xreg = St_reg )
127
18
save.Object ( ARMA00 ,"/ Volumes / USB / Parameterestimater / Storage_Reg.txt " )
19
# Benytter kun residulaer fra 29 . marts. 2010
20
reskort ← resid ( ARMA00 ) [133: NROW ( resid ( ARMA00 ))]
21 22
AIC ← matrix ( numeric () , nrow =5 , ncol =5)
23
BIC ← matrix ( numeric () , nrow =5 , ncol =5)
24
rownames ( AIC ) ← paste (" AR " ,1:5 , sep = "")
25
colnames ( AIC ) ← paste (" MA " ,1:5 , sep = "")
26
rownames ( BIC ) ← paste (" AR " ,1:5 , sep = "")
27
colnames ( BIC ) ← paste (" MA " ,1:5 , sep = "")
28
for (j in 1:5) {
29
for ( i in 1:5) { fit ← Arima ( reskort ,c(i ,0 , j) , include.mean = FALSE )
30 31
AIC [i ,j]
←
fit $ aic
32
BIC [i ,j]
←
fit $ bic
33
}
34
}
35 36
NY ← Arima ( reskort ,c (2 ,0 ,1) , include.mean = FALSE )
37
save.Object (NY , "/ Volumes / USB / Parameterestimater / Storage_ARMA_Kort.txt ")
38
summary ( NY )
39 40
PX ← as.numeric ( predict ( ARMA00 , newxreg = pred_cos_sin [, -5 ]) $ pred )
Listing 3: Lagermodel
Oliemodel 1
dat ← read.csv ( "/ Volumes / USB / April 2014 Data / Oil_to_R April2014.csv " , header = TRUE , sep = ";" , dec =" ,")
128
2
yt ← as.matrix ( dat [ ,2:9])
3
ttm ←
as.matrix ( dat [ ,17:24])
4
dt
dat [ ,32]
←
5 6
adf.test ( yt [ ,1])
7
# stationaer
8
CO1_afkast ← diff ( log ( yt )) [ ,1]
9 10
adf.test ( CO1_afkast ); acf ( ts ( CO1_afkast )) ; pacf ( ts ( CO1_afkast ))
11 12
# ## Fitter MA1
13
MA1 ← Arima ( CO1_afkast , order =c (0 ,0 ,1) ) ; summary ( MA1 )
14 15
res ← residuals ( MA1 )
16 17
# ## Finder min BIC
18
bics ← matrix ( nr =5 , nc =5)
19
colnames ( bics ) ← rownames ( bics ) ← 0:4
20
for (i in 2:5) for (j in 1:5) {
21
spec ← as.formula ( paste ( '∼garch ( ',i-1 , ',' ,j-1 , ') '))
22
bics [i ,j ] ← garchFit ( spec , res , trace =F ) @fit $ ics [2]
23
}
24
which ( bics == min ( bics , na.rm = TRUE ))
25 26
# ## Finder min BIC
27
bics ← matrix ( nr =5 , nc =5)
28
colnames ( bics ) ← rownames ( bics ) ← 0:4
29
for (i in 2:5) for (j in 1:5) {
30
spec ← as.formula ( paste ( '∼garch ( ',i-1 , ',' ,j-1 , ') '))
31
bics [i ,j ] ← garchFit ( spec , res , trace =F , cond.dist =" sstd ") @fit $ ics
129
[2] 32
}
33 34
which ( bics == min ( bics , na.rm = TRUE ) )
35 36
# ## Fitter bedste Garch
37
g← garchFit (∼garch (1 ,1) ,res , trace =F , include.mean = FALSE )
38
summary ( g)
39 40
# ## Residualer ikke normalfordelt
41
par ( mfrow =c (1 ,2) )
42
res_g ← residuals (g , standardize =T )
43 44
# Are they normally distributed ? No way.
45
qqnorm ( res_g , xlim =c( -6 ,4) ); qqline ( res_g )
46 47
# ## Finder maksimum likelihood estimat for t-fordelt data
48
ll ← function ( scale , df ) -sum ( dt ( residuals (g ,T) / scale , df =df , log =T) -log ( scale ) )
49
maxLike ← mle (ll , start = list ( scale = sd ( residuals (g , T)) ,df =5) )
50
degf ← maxLike@coef [2]
51 52
# ## Dette ser bedre ud
53
plot ( qt ( ppoints ( residuals (g ,T) ) ,df = degf ) , sort ( residuals (g , T)) , xlab = ' Theoretical Quantiles ', ylab = ' Sample Quantiles ', main = 't Q-Q Plot ')
54
my.qt ← function ( x) qt (x , df = degf )
55
qqline ( residuals (g ,T ) , distribution = my.qt )
56 57
# ## Fitter denne og ser , at alt er signifikant
58
g.t ← garchFit (∼garch (1 ,1) ,res , trace =F , cond.dist = ' std ')
130
59
summary ( g.t )
60 61
# ################# Skewed t
62
g.t.s ← garchFit (∼garch (1 ,1) ,res , trace =F , cond.dist = ' sstd ', includ.skew = TRUE )
63
shape ← g.t.s@fit $ matco [ " shape " ,1]
64
skew ← g.t.s@fit $ matco [" skew " ,1]
65 66
# ## Dette ser bedre ud
67
plot ( qsstd ( ppoints ( residuals ( g.t.s ,T )) , nu = shape , xi = skew ) , sort ( residuals ( g.t.s , T)) , xlab = ' Theoretical Quantiles ', ylab = ' Sample Quantiles ', main = 't skewed Q-Q Plot ', xlim =c( -6 ,6) )
68
my.qt ← function ( x) qsstd (x , nu = shape , xi = skew )
69
qqline ( residuals ( g.t.s ,T ) , distribution = my.qt )
70 71
# ## Fitter det hele paa en gang med skeewed t fordelt
72
summary ( arma.garch ← garchFit (∼arma (0 ,1) + garch (1 ,1) , CO1_afkast , trace =F , cond.dist = ' sstd ') )
73
res.a.g ← residuals ( arma.garch ,T )
74 75
# ############ Simulation
76
est ← arma.garch@fit $ matco [ ,1]
77 78
spec_sim = garchSpec ( model = list ( mu = est [" mu "] , ma = est [" ma1 "], omega = est [" omega " ], alpha = est [ " alpha1 "], beta = est [ " beta1 "], skew = est [" skew "] , shape = est [" shape " ]) , cond.dist = " sstd " )
79
garchSim ( spec_sim , n =365)
80
plot.ts ( garchSim ( spec_sim , n =365) )
81 82
oilsim ← function () {
131
83
g_sim ← garchSim ( spec_sim , n =365)
84
oil_sim ← yt [ NROW ( yt ) ,1] * exp ( -cumsum ( g_sim ) )
85
return ( oil_sim )
86
}
87 88
# 180 dages gennemsnit over de 180 sidste dage
89
# Med 130 arbejdsdage
90
moving_average ← function () {
91
frontM ← c( yt [( NROW ( yt ) -128 ): NROW ( yt ) ,1] , oilsim () )
92
M ← numeric (365)
93
for (k in 0:364) M[k +1] ← mean ( frontM [( k +1) :( k +130) ])
94
return (M )
95
}
Listing 4: Oliemodel
Simulation af fejlled 1
fejl_sim ← function ( initial = tail ( resid ( model1 ) ,5) ){
2
x ← numeric (365)
3
noise ← rghyp (365 , NIG )
4
x [1] ← AR5 $ coef %*% initial [5:1] + noise [1]
5
x [2] ← AR5 $ coef %*% c(x [1] , initial [5:2]) + noise [2]
6
x [3] ← AR5 $ coef %*% c(x [2] , x [1] , initial [5:3]) + noise [3]
7
x [4] ← AR5 $ coef %*% c(x [3] , x [2] , x [1] , initial [5:4]) + noise [4]
8
x [5] ← AR5 $ coef %*% c(x [4] , x [3] , x [2] , x [1] , initial [5]) + noise [5]
9 10
for (i in 6:365)
11
return (x )
12
x[i ] ← AR5 $ coef %*% x [( i-1 ) :( i-5 )] + noise [i]
}
Listing 5: Simulation af fejlled
132
Simulation af olie 1
oilsim ← function ( slutolie = yt [ NROW ( yt ) ]) {
2
g_sim ← garchSim ( spec_sim , n =365)
3
oil_sim ←
4
return ( oil_sim )
5
}
slutolie * exp ( -cumsum ( g_sim ) )
Listing 6: Simulation af olie 1
moving_average ← function ( oliesidst = yt [ NROW ( yt )] , index = NROW ( yt ) ){
2
frontM ← c( yt [( index-128 ): index ], oilsim ( oliesidst ))
3
M ← numeric (365)
4
for (k in 0:364) M[k +1] ← mean ( frontM [( k +1) :( k +130) ])
5
return (M )
6
}
Listing 7: Simulation af moving average
Simulation af Lager 1
sim_gas_storage ← function ( sidste_storage_pct , new_reg , storage_gns , Future = TRUE ){
2 3
# Modellen er foretaget paa logit til procent.
4
PRE ← as.numeric ( predict ( Storage.ARMAreg , newxreg = new_reg )$ pred )
5
sim1 ← 1/ (1+ exp (-( simulate.Arima ( Storage.ARMAkort ,53 , future = Future )+ PRE )))
6
# Lineaer interpolerer mellem de simulerede mandage. int_sim ← approx (x =c (0 , seq (7 ,371 ,7) ) ,y=c ( sidste_storage_pct , sim1
7
) ,n =372) $y [2:366] int_sim-storage_gns
8 9
}
Listing 8: Simulation af lager 133
Simulation af Temperatur 1
Simulation_Temp ← function ( future = TRUE ) { resultat ← simulate.Arima ( Temp_model , nsim = NROW ( ny_reg ) , xreg =
2
ny_reg , future = future ) return ( resultat )
3 4
}
Listing 9: Simulation af temperatur 1
Simulation_W_shocks ← function ( gns_W_shocks , Future = TRUE ) {
2
resultat ← gns_W_shocks - simulate.Arima ( Temp_model , nsim = NROW ( ny_reg ) , xreg = ny_reg , future = Future )
3
return ( resultat )
4
}
Listing 10: Simulation af temperaturchok
Samlet simulation 1
sti ← "/ Volumes / USB /"
2
sti ← "F :/"
3
nem ← function ( path_end ) { paste ( sti , path_end , sep =" ")}
4 5
# ----------------- OIL
6
source ( nem ( ' Simulation / oilsim.R ') )
7
source ( nem ( ' Simulation / moving_average.R '))
8
yt ← read.csv ( nem (" April 2014 Data / Oil_to_R April2014.csv ") , header = TRUE , sep = ";" , dec =" ,") [ ,2]
9
est ← loadObject ( nem ( " Parameterestimater / Oil_Garch.txt ") )
10
spec_sim = garchSpec ( model = est , cond.dist = " sstd " )
11
OlieMovDat ← read.csv ( nem (" April 2014 Data / Samlet data APRIL 2014 .csv " ) , header = TRUE , sep =" ;" , dec =" ,")[ ,c (2 ,16) ]
134
12
OlieDat ← read.csv ( nem (" April 2014 Data / Samlet data APRIL 2014 .csv " ) , header = TRUE , sep ="; " , dec =" ,") [,c (2 ,24) ]
13 14
# ------------------- Temperatur
15
source ( nem ( ' Simulation / Simulation_W_shocks.R '))
16
source ( nem ( ' Simulation / Simulation_Temp.R ') )
17
ny_reg ← read.table ( nem ( " Simulation / Ny_reg_temp.txt " ) , header = TRUE )
18
Temp_model ← loadObject ( nem ( " Parameterestimater / Temp_AR4.txt "))
19
Temp_sim_input ← read.csv ( nem (" April 2014 Data / sim_temp_til_r.csv ") , header = TRUE , sep =" ;" , dec =" ,")
20
W_shocks1 ← read.csv ( nem (" April 2014 Data / Samlet data APRIL 2014 .csv " ) , header = TRUE , sep =" ;" , dec =" ,")[ ,c (2 ,17) ]
21
TempDat ← read.csv ( nem (" April 2014 Data / Temperatur-model til R.csv " ) , header = TRUE , sep ="; " , dec =" ,") [ ,1]
22
TempDat1 ← TempDat [( NROW ( TempDat ) -NROW ( W_shocks1 [ ,2]) +1) : NROW ( TempDat ) ]
23 24
# ------------------- Deterministisk data gas-year-2015
25
Det_dat ← read.csv ( nem (" Simulation / Deterministik_data_til_sim_R.csv ") , header = TRUE , sep =" ;" , dec = " ," )
26 27
# -------------------- Gasstorage
28
new_reg
← tail ( read.csv ( nem (" April 2014 Data / Storage
sin-cos-reg.csv ") , header = TRUE , sep = ";" , dec =" ,") [, -c (1:2 ,7:12) ] ,371) [ seq (7 ,371 ,7) ,-5 ] 29
Storage.ARMAreg ← loadObject ( nem ( " Parameterestimater / Storage_Reg.txt "))
30
Storage.ARMAkort ← loadObject ( nem ( " Parameterestimater / Storage_ARMA_Kort.txt "))
31
135
32 33
sidste_storage_pct ← read.csv ( nem (" April 2014 Data / Storage sin-cos-reg.csv ") , header = TRUE , sep = ";" , dec =" ,") [2388 ,12]
34
storage_gns
← read.csv ( nem (" Simulation / Storage_GNS_input_sim.csv ")
, header = TRUE , sep =" ;" , dec = " ," ) [ ,2] 35
source ( nem ( ' Simulation / sim_gas_storage.R ') )
36 37
sim_gas_storage ( sidste_storage_pct , new_reg , storage_gns )
38
plot.ts ( sim_gas_storage ( sidste_storage_pct , new_reg , storage_gns ))
39 40
# -------------------- Load parameter estimater og data til Gas model
41
source ( nem ( ' Simulation / fejl_sim.R '))
42
model1 ← loadObject ( nem ( " Parameterestimater / Model1.txt ") ); summary ( model1 )
43
dat1 ← read.csv ( nem (" April 2014 Data / Samlet data APRIL 2014 .csv ") , header = TRUE , sep =" ;" , dec =" ," ) [ ,2:3]
44
AR5 ← loadObject ( nem ( " Parameterestimater / Gas_AR5.txt " ))
45
NIG ← loadObject ( nem ( " Parameterestimater / NIG.txt " ))
46
new_reg1
← read.csv ( nem (" April 2014 Data / Storage
sin-cos-reg.csv ") , header = TRUE , sep = ";" , dec =" ,") [2389:2753 , -c (1:2 ,7:13) ] 47 48
# -------------------- Dan simulationsdata
49
Storage_ny_hat
← read.csv ( nem (" Simulation /
Storage_GNS_input_sim.csv " ) , header = TRUE , sep =" ;" , dec =" ,") [ ,2] 50
nydata ← function () {
51
NYDAT ← cbind ( Det_dat , moving_average () , as.vector ( Simulation_W_shocks ( Temp_sim_input [ ,3]) ) , new_reg1 , Storage_ny_hat , sim_gas_storage ( sidste_storage_pct , new_reg , storage_gns ))
52
names ( NYDAT )[ c (6:7 ,13) ] ← c( " X.180.dages.MA..med.interpolation " ,"
136
W_shocks " ," Storage_andel_ugentlig " ) 53
return ( NYDAT )
54
}
55 56
# -------------------- GAS simulation og plot
57
exp ( predict ( model1 , newdata = nydata () ) + fejl_sim () )
58 59
# -------------------- GAS simulation data til simpelt eksempel
60
x ← nydata ()
61
Y← cbind (x , exp ( predict ( model1 , newdata = x) + fejl_sim () ) )[ seq (1 ,365 ,60) ,-c (1:5 ,8:12) ]
62
SIM ← array (, dim =c (8 ,7 ,4) , dimnames =c ( ,) )
63 64
for (b in 1:8) {
65
x ← nydata ()
66
Y← cbind (x , exp ( predict ( model1 , newdata = x) + fejl_sim () ) )[ seq (1 ,365 ,60) ,-c (1:5 ,8:12) ]
67
Y ← t( Y)
68
for (j in 1:4) {
69
SIM [b ,,j ] ← Y[j ,]
70
}}
71 72
rownames ( SIM ) ← c( " Oliepris " ," Temperaturshock " ," Lagerniveau " ," Lagershock " ," Gaspris ")
73
SIM ← list ( Oliepris = SIM [ , ,1] , Temperaturshock = SIM [ , ,2] , Lagershock = SIM [ , ,3] , Gaspris = SIM [ , ,4])
74
# saveObject ( SIM , nem (" Simulation / SIMeksempel.txt ") )
Listing 11: Samlet simulation
Udglatning af forwardkurven 137
1
C ← read.csv ( "/ Volumes / USB / Data / Grunddata / GASfuturesKalibrering.csv " , sep ="; " , dec =" ,") [ ,2:9]
2
n ← 8
3
m ← 8
4
H ← matrix (0 , nrow =5* n , ncol =5* n)
5
A ← matrix (0 , nrow =3*( n-1 )+m +1+1 , ncol =5* n )
6
t ← c (0 ,31 ,62 ,92 ,123 ,154 ,184 ,276 ,365)
7 8
# Note quarter month problem
9
pris ← C [8:1]
10
pris1 ← pris
11
pris1 [6] ← 3*( pris [6] - ( sum ( pris [4:5]) )/ 3)
12
names ( pris1 ) [6] ← " sep.14 "
13
D ← cbind (t ( pris1 ) ,t [-( n +1) ], t[ -1 ], diff (t) )
14 15
for (j in 1: n ){
16
h ← matrix (c (144 / 5* diff ( t ∧ 5) [ j ] ,18* diff (t ∧ 4) [ j ] ,8* diff (t ∧ 3) [ j ] ,0 ,0 ,18* diff ( t ∧ 4) [ j ] ,12* diff (t ∧ 3) [ j ] ,6* diff (t ∧ 2) [ j ] ,0 ,0 ,8* diff (t ∧ 3) [ j ] ,6* diff (t ∧ 2) [ j ] ,4* diff (t) [j] , rep (0 ,12) ) , nrow =5)
17
H [(5*( j-1 ) +1) :(5*( j-1 ) +5) ,(5*( j-1 ) +1) :(5*( j-1 ) +5) ] ← h
18
}
19 20
for (j in 1:( n-1 ) ){
21
startC ← ( j-1 ) *5+1
22
startR ← ( j-1 ) *3+1
23
A[ startR , startC :( startC +9) ] ← c( -t [ j +1] ∧ (4:0) , t[j +1] ∧ (4:0) )
24
A[ startR +1 , startC :( startC +9) ] ← c( - (4:0) *t[ j +1] ∧ (c (3:0 ,0) ) ,(4:0) * t[ j +1] ∧ (c (3:0 ,0) ) )
25
A[ startR +2 , startC :( startC +9) ] ← c( -(c (12 ,6 ,2 ,0 ,0) )*t [j +1] ∧ (c (2:0 ,0 ,0) ) ,(c (12 ,6 ,2 ,0 ,0) ) *t[ j +1] ∧ (c (2:0 ,0 ,0) ) )
138
26
}
27 28
A [( n-1 ) *3+1 ,(( n-1 ) *5+1) :(( n-1 ) *5+5) ] ← (4:0) *t [n] ∧ (c (3:0 ,0) )
29
A [3*( n-1 )+m +1+1 ,5] ← 1
30 31 32
for (k in 1: m ){
33
for (j in 0:4) {
34
A [( n-1 ) *3+1+ k ,( j +1+5*( k-1 )) ] ← (1 / (5 -j )) * diff ( t ∧ (5 -j )) [k] / diff ( t)[ k]
35
print ( cbind (k ,j ,(1 / (5 -j )) , diff (t ∧ (5 -j )) [k] , diff ( t)[ k ]) )
36
}
37
}
38
Giant ← rbind ( cbind (2* H , t(A )) , cbind (A , matrix (0 , ncol =31 , nrow =31) ) )
39
b ← rep (0 ,8* n+ m-2 +1)
40
b [(8* n-1 ) :(8* n+ m-2 )] ← D [ ,1]
41
b [8* n + m-2 +1] ← 21 .40
42
coef ← qr.solve ( Giant ,b , tol =1 E-115 )
43
ti ← 1:365
44
funk ← function ( ti , coef ,m ){
45
res ← numeric ( length ( ti ))
46
for (i in 1: length ( ti )) {
47
for ( j in 1: m ){
48
if ( t[j ]
49
}
50
}
51
}
52
return ( res )
139
53
}
Listing 12: Forward smooth
Kalibrering til forwardkurven 1
sti ← "F :/"
2
nem ← function ( path_end ) { paste ( sti , path_end , sep =" ")}
3
P ← read.csv ( nem (" Simulation / Simulerede_stier / Gas_uden_outliers_uden_mod.csv " ) , sep = ";" , dec =" ," , header = FALSE )
4
sf ← as.matrix ( read.table ( nem ( " Simulation / Simulerede_stier / smoothforward.txt ")))
5
P ← t( P)
6
Pmod ← P *0
7
rsp ← rowSums ( P)
8
for (t in 1:365) { Pmod [t ,] ← ( NCOL (P) * P[t ,] / rsp [t ]) * sf [t ,1]
9 10
}
Listing 13: Kalibrering til forwardkurven
13.5 Kode til at danne resultater LSMC Gaslagerkontrakt 1
sti ← "F :/"
2
nem ← function ( path_end ) { paste ( sti , path_end , sep =" ")}
3 4
# ######## Input
5
P ← read.table ( nem (" Simulation / Simulerede_stier / GASSIMmodificeretstor.txt " ) , sep = ";" )
6
#P ← loadObject ( nem ( " Simulation / SIMeksempel.txt " ))$ Gaspris [, -c (2 ,4 ,6) ]
7
140
8
d ← 1 .02 ∧ ( -1 * 60 / 360)
9 10
(i ← 1 .65 )
11
(w ← 2 .2 )
12 13
( Vmax ← 3)
14
( Vmin ← 0)
15
( precision ← 2)
16
( ai ← 0)
17
( bi ← 0)
18
( aw ← 0)
19
( bw ← 0)
20 21
# floor function
22
nyfloor ← function ( x){
23
if ( abs ( x-floor ( x)) >1 -10 ∧ ( -9 ) ) { res ← floor (x ) + 1}
24
if ( abs ( x-floor ( x))≤1 -10 ∧ ( -9 )) { res ← floor (x ) }
25
return ( res )
26
}
27 28
nyceiling ← function ( x){
29
if ( abs ( x-ceiling (x )) >1 -10 ∧ ( -9 ) ) { res ← ceiling (x) - 1}
30
if ( abs ( x-ceiling (x ))≤1 -10 ∧ ( -9 )) { res ← ceiling ( x) }
31
return ( res )
32
}
33 34
# #########
35
( Vslut ← Vmin )
36
( dv ← w/ precision )
37
v ← seq ( Vmin , Vmax , dv ) ; if (v[ length (v) ]6=Vmax ) { v ← c(v , Vmax ) }; v
141
38
Nsim ← NROW ( P)
39
Ntid ← NCOL ( P)
40
Nm
← length (v )
41 42
# ## Grid
43
( th ← nyfloor (( Vmax-Vslut ) /w) )
44
( rand ← c( rep ( Vmax , Ntid-th ) ,w *( th :1) ) )
45 46
# Tvungne widrawal tidspunkter
47
(j ← seq (( ( Ntid +1) - th ) , Ntid ,1) )
48 49
# payoff matricer for forskellige volume niveauer
50
( VM ← array (, dim =c( Nsim , Ntid , Nm ) ) )
51
# Continuation value
52
(C ← array (, dim =c( Nsim , Ntid , Nm ) ) )
53 54
# Vmin matricen
55
( VM [, Ntid ,1] ← 0)
56 57
# Sidste vaerdi
58
for (k in 1: Nm ){ VM [ , Ntid ,k ] ← min (v [k] ,w) * ( P[, Ntid ] *(1 -aw ) - bw )}
59
# Vaerdi over rand
60
for (t in 1:( Ntid-1 ) ){
61
for (k in 1: Nm ){ if ( v[k ]> rand [t ]) { VM [,t ,k]← w * ( P [,t ] *(1 -aw ) - bw ) }
62 63
}}
64 65
# Tvungne payoffs
66
lj ← length (j )
142
67
for ( id in 1: lj ){ VM [ ,j[ id : lj ] , precision *( lj-length ( j [1: id ]) +1) +1] ← w * ( P [,j [ id
68
: lj ]] *(1 -aw ) - bw ) 69
}
70
VMtemp ← VM
71
Y ← matrix (, ncol =Nm , nrow = Nsim )
72 73
# Step 1 fyld matricerne , saa lige mange soeler har tal.
74
for (t in ( Ntid-1 ) :1 ){
75
ptm ← proc.time ()
76
for (k in 1: Nm ){Y [,k ] ←
as.matrix ( VM [ ,( t +1) : Ntid , k ]) %*% d ∧ (1:(
Ntid-t ) ) } 77 78
C ← fitted ( lm ( Y
∼ P[ ,t ]+ I( P[, t] ∧ 2) ) )
79 80
# k = 1
81
k = 1
82
imsd ← k+ nyceiling ( min (i , rand [t +1] -v [ k ]) / dv ) -1
# linds op og saa
ned ide 83
imsu ← imsd + 1
84
xi ← (v [k] + min (i , rand [ t +1] -v [k ]) [ imsd ])
85
Decisionmatrix ← cbind (C [,k ], *(1 + ai ) + bi
86
;
) +
- v [ imsd ] ) /
min (dv ,v[ imsu ] -v
print ( xi ) - min (i , rand [ t +1] -v [ k ]) * ( P[, t]
xi *C[ , imsu ] + (1 -xi )*C [, imsd ] )
l ← apply ( Decisionmatrix ,1 , which.max ) -1
87 88
# inject payoff or do nothing
89
VM [ ,t , k] ← -min (i , rand [t +1] -v [k ]) *( P[, t] *(1 + ai ) + bi
90
for (b in 1: Nsim ) {
91
VMtemp [b ,( t +1) : Ntid , k] ←
) *l
xi *( l[ b ]==1) * VM [b ,( t +1) : Ntid ,k + l[ b ]*(
143
l[ b ]==1) *( imsu-k )]
+
(1 -xi ) *( l [b ]==1) * VM [b ,( t +1) : Ntid , k + l [b ]*( l [b
92
]==1) *( imsd-k )]
+
(l [b ]==0) * VM [b ,( t +1) : Ntid , k ]
93 94 95
}
96
# ######### k under rand
97
for ( k in 2:( Nm-1 )) {
98
# Ved Aendring a Grid skal ses paa C [, matrice skift ] , som ikke noedvendigvis er 1.
99
if ( v[k ]≤rand [ t +1]) {
100
# i : injection , w : widrawal , ms : matriceskift , u:up , d : down
101
imsd ← k+ nyceiling ( min (i , rand [t +1] -v [ k ]) / dv ) - 1
102
imsu ← imsd + 1
103
wmsd ← k-nyceiling ( min (w , v[k ] -Vmin )/ dv )
104
wmsu ← wmsd + 1
105 106
# up weigths
107
xi ← (v [k] + min (i , rand [ t +1] -v [k ]) [ imsd ])
108
min (dv ,v[ imsu ] -v
; # print ( xi )
xw ← (v [k] - min (w ,v[ k] -Vmin ) ])
- v [ imsd ] ) /
- v[ wmsd ] ) /
min (dv ,v [ wmsu ] -v [ wmsd
; # print ( xw )
109 110
Decisionmatrix ← cbind ( min (w ,v [k] -Vmin ) * ( P[,t ] *(1 - aw ) - bw ) + xw * C[, wmsu ] + (1 -xw ) *C[ , wmsd ] , C[ ,k] , ( P [,t ] *(1 + ai ) + bi
) + xi * C[, imsu ]
- min (i , rand [t +1] -v [k ]) * + (1 -xi )*C[ , imsd ] )
111 112
l ← apply ( Decisionmatrix ,1 , which.max ) - 2
113 114
VM [ l == -1 ,t ,k ]
←
min (w ,v[ k] -Vmin ) * ( P [l == -1 , t] *(1 - aw ) - bw )
144
115
VM [ l ==0 ,t ,k]
←
0
116
VM [ l ==1 ,t ,k]
←
- min (i , rand [t +1] -v [k ]) * ( P[ l ==1 , t ] *(1 + ai )
+ bi
)
117
for (b in 1: Nsim ) {
118
VMtemp [b ,( t +1) : Ntid , k] ← ( xw *( l[b ]== -1 )+ xi *( l[ b ]==1) ) * VM [b ,( t +1) : Ntid ,k + l [b ]*( l[ b ]== -1 ) *( k-wmsu ) + l [b ]*( l [b ]==1) *( imsu-k ) ]
+
((1 -xw ) *( l[ b ]== -1 ) +(1 -xi ) *( l[ b ]==1) ) * VM [b ,( t +1)
119
: Ntid ,k + l[b ]*( l[b ]== -1 ) *( k-wmsd ) + l[ b ]*( l[ b ]==1) *( imsd-k )]
+
(l [b ]==0) * VM [b ,( t +1) : Ntid , k ]
120 121
}}
122
if ( rand [t +1]
123
wmsd ← k-nyceiling ( min (w , v[k ] -Vmin )/ dv )
124
wmsu ← wmsd + 1
125
xw ← (v [k] - min (w ,v[ k] -Vmin ) ])
;#
- v[ wmsd ] ) /
min (dv ,v [ wmsu ] -v [ wmsd
print ( xw )
126 127
Decisionmatrix ← cbind ( min (w ,v [k] -Vmin ) * ( P[,t ] *(1 - aw ) - bw ) + xw * C[, wmsu ] + (1 -xw ) *C[ , wmsd ] ,
(v [k] -rand [ t +1]) * ( P [,t] *(1 -
aw ) - bw ) + C[, which ( v == rand [t +1]) ] ) 128 129
l ← apply ( Decisionmatrix ,1 , which.max ) - 2
130 131
VM [ l == -1 ,t ,k ]
←
min (w ,v[ k] -Vmin ) * ( P [l == -1 , t] *(1 - aw ) - bw )
132
VM [ l ==0 ,t ,k]
←
133
for (b in 1: Nsim ) {
134
VMtemp [b ,( t +1) : Ntid , k] ← VM [b ,( t +1) : Ntid ,k+
(v [k] -rand [ t +1]) * ( P[ l ==0 , t ] *(1 - aw ) - bw )
l[ b ]*( l[ b ]== -1 ) * min (w ,v[ k] -Vmin )/ dv - (l[ b ]==0)
135
*( v [k] -rand [ t +1]) / dv ] 136
}}}
145
137
# ######### Sidste k
138
if ( v[ Nm ]== rand [ t ]) {
139
k= Nm
140 141
wmsd ← k - 1 - nyceiling (( w- (v[ k] -v [ k-1 ]) )/ dv
142
wmsu ← wmsd + 1
143
xw ← (v [k] - min (w ,v[ k] -Vmin )
- v[ wmsd ] ) /
)
dv
144 145
Decisionmatrix ← cbind ( min (w ,v [k] -Vmin ) * ( P[,t ] *(1 - aw ) - bw ) + xw * C[, wmsu ] + (1 -xw ) *C[ , wmsd ] , (v [k] -rand [ t +1]) * ( P[ ,t] *(1 aw ) - bw ) + C[, which ( v == rand [t +1]) ])
146 147
l ← apply ( Decisionmatrix ,1 , which.max ) - 2
148 149
VM [ l == -1 ,t ,k ]
←
min (w ,v[ k] -Vmin ) * ( P [l == -1 , t] *(1 - aw ) - bw )
150
VM [ l ==0 ,t ,k]
←
151
for (b in 1: Nsim ) {
152
VMtemp [b ,( t +1) : Ntid , k] ← xw * VM [b ,( t +1) : Ntid ,k+
(v [k] -rand [ t +1]) * ( P[ l ==0 , t ] *(1 - aw ) - bw )
l[ b ]*( l[ b ]== -1 ) *( k-wmsu ) - ( l[b ]==0) *( k-which (v
153
== rand [t +1]) )] + (1 -xw ) * VM [b ,( t +1) : Ntid ,k + l[ b ]*( l[ b ]== -1 ) *( k-wmsd ) - ( l[b ]==0) *( k-which (v
154
== rand [t +1]) )] 155
}}
156
# Updating step
157
VM [ ,( t +1) : Ntid ,]← VMtemp [ ,( t +1) : Ntid ,]
158
Timemeasure [t]← proc.time () [ " elapsed "] -ptm [" elapsed "]
159
}
Listing 14: LSMC Gaslager
146
LSMC Hedge 2014-2015 1
nem ← function ( path_end ) { paste ( sti , path_end , sep =" ")}
2 3
# ######## Input
4
ptm ← proc.time ()
5
Storage ← read.csv ( nem (" Simulation / Simulerede_stier / Storage.csv ") , sep =" ;" , dec =" ," , header = FALSE )
6
Temp ← read.csv ( nem (" Simulation / Simulerede_stier / Temperatur.csv ") , sep =" ;" , dec =" ," , header = FALSE )
7
Oil ← read.csv ( nem (" Data / Simulationsdata / Olie30000.csv " ) , sep =" ;" , dec =" ," , header = FALSE )
8
P1 ← read.csv ( nem (" Simulation / Simulerede_stier / Gas.csv " ) , sep = ";" , dec =" ," , header = FALSE )
9
Dat15 ← read.csv ( nem (" April 2015 Data / Behandlet data / Samlet_data_14-15.csv ") , sep ="; " , dec = " ," , header = TRUE )
10
proc.time () [" elapsed "] -ptm [" elapsed "]
11 12
ptm ← proc.time ()
13
Antalsim ← 100
14
Indikator ← c( rep ( FALSE ,91) , rep ( TRUE ,184) , rep ( FALSE ,90) )
15
P← -P1 [1: Antalsim , Indikator ] ; NROW (P ); NCOL ( P)
16
Olie ← Oil [1: Antalsim , Indikator ]; NROW ( Olie ) ; NCOL ( Olie )
17
Temperatur ← Temp [1: Antalsim , Indikator ]; NROW ( Temperatur ) ; NCOL ( Temperatur )
18
Lager ← Storage [1: Antalsim , Indikator ]; NROW ( Lager ); NCOL ( Lager )
19
proc.time () [" elapsed "] -ptm [" elapsed "]
20 21
P ← rbind (P , -t ( Dat15 $ Gas )[ Indikator ])
22
Olie ← rbind ( Olie ,t( Dat15 $ X180.dages.MA.med.interpolation )[ Indikator ])
147
23
Temperatur ← rbind ( Temperatur ,t( Dat15 $ W_shocks )[ Indikator ])
24
Lager ← rbind ( Lager ,t ( Dat15 $ Storage_andel_ugentlig )[ Indikator ])
25 26
# ### ALLE STATE VARIABLE
27
d ← 1 .0 ∧ ( -1 * 60 / 360)
28
(i ← 3000)
29
(w ← 1500)
30
( Vmax ← 100000)
31
( Vmin ← 0)
32
( precision ← 1)
33
( ai ← 0)
34
( bi ← 0)
35
( aw ← 0)
36
( bw ← 0)
37 38
# floor function
39
nyfloor ← function ( x){
40
if ( abs ( x-floor ( x)) >1 -10 ∧ ( -9 ) ) { res ← floor (x ) + 1}
41
if ( abs ( x-floor ( x))≤1 -10 ∧ ( -9 )) { res ← floor (x ) }
42
return ( res )
43
}
44 45
nyceiling ← function ( x){
46
if ( abs ( x-ceiling (x )) >1 -10 ∧ ( -9 ) ) { res ← ceiling (x) - 1}
47
if ( abs ( x-ceiling (x ))≤1 -10 ∧ ( -9 )) { res ← ceiling ( x) }
48
return ( res )
49
}
50 51
# #########
52
( Vslut ← Vmin )
148
53
( dv ← w/ precision )
54
v ← seq ( Vmin , Vmax , dv ) ; if (v[ length (v) ]6=Vmax ) { v ← c(v , Vmax ) }; v
55
Nsim ← NROW ( P)
56
Ntid ← NCOL ( P)
57
Nm
← length (v )
58 59
# ## Grid
60
( th ← nyfloor (( Vmax-Vslut ) /w) )
61
( rand ← c( rep ( Vmax , Ntid-th ) ,w *( th :1) ) )
62 63
# Tvungne widrawal tidspunkter
64
(j ← seq (( ( Ntid +1) - th ) , Ntid ,1) )
65 66
# payoff matricer for forskellige volume niveauer
67
( VM ← array (, dim =c( Nsim ,2 , Nm ) ) )
68
# Beta array
69
( beta ← array (, dim =c (5 , Ntid-1 , Nm ) ))
70
# Vmin matricen
71
( VM [ ,2 ,1] ← 0)
72 73
VMtemp ← array (, dim =c( Nsim ,1 , Nm ) )
74
Y ← matrix (, ncol =Nm , nrow = Nsim )
75
Timemeasure ← numeric ( Ntid-1 )
76 77
# Step 1 fyld matricerne , s å lige mange sø jler har tal.
78
for (t in ( Ntid-1 ) :1 ){
79
ptm ← proc.time ()
80
for (k in 1: Nm ){ Y[,k ] ←
81
X ← as.matrix ( cbind (1 , P [,t ], Olie [,t ], Lager [, t], Temperatur [,t ]) )
82
fit ← lm.fit (X [ -NROW (P) ,] , Y[ -NROW (P ) ,] )
VM [ ,2 , k ]* d }
149
83
C ← X %*% fit $ coefficients
84
beta [ ,t ,] ←
fit $ coefficients
85 86
# k = 1
87
k = 1
88
imsd ← k+ nyceiling ( min (i , rand [t +1] -v [ k ]) / dv ) -1
# linds op og s å
c ned idà 89
imsu ← imsd + 1
90
xi ← (v [k] + min (i , rand [ t +1] -v [k ]) [ imsd ])
91
;
Decisionmatrix ← cbind (C [,k ], *(1 + ai ) + bi
92
) +
- v [ imsd ] ) /
min (dv ,v[ imsu ] -v
# print ( xi ) - min (i , rand [ t +1] -v [ k ]) * ( P[, t]
xi *C[ , imsu ] + (1 -xi )*C [, imsd ] )
l ← apply ( Decisionmatrix ,1 , which.max ) -1
93 94
# inject payoff or do nothing
95
VM [ ,1 , k] ← -min (i , rand [ t +1] -v [k ]) *( P[, t] *(1 + ai ) + bi
96
l0 ← l ==0
97
l1 ← l ==1
) *l
98 99 100
for (b in 1: Nsim ) { VMtemp [b ,1 , k] ←
xi *( l1 [b ]) * VM [b ,2 , k + l[b ]*( l1 [b ]) *( imsu-k )] (1 -xi ) *( l [b ]==1) * VM [b ,2 , k + l[b ]*( l1 [b ]) *(
101
imsd-k )]
+
( l0 [b ]) * VM [b ,2 , k ]
102 103
}
104
# ######### k under rand
105
for ( k in 2:( Nm-1 )) {
106
# Ved æ ndring a Grid skal ses på C[, matrice skift ], som ikke nø dvendigvis er 1.
107
if ( v[k ]≤rand [ t +1]) {
150
+
108
# i : injection , w : widrawal , ms : matriceskift , u:up , d : down
109
imsd ← k+ nyceiling ( min (i , rand [t +1] -v [ k ]) / dv ) - 1
# linds op og så
c ned idà 110
imsu ← imsd + 1
111
wmsd ← k-nyceiling ( min (w , v[k ] -Vmin )/ dv )
112
wmsu ← wmsd + 1
113 114
# up weigths
115
xi ← (v [k] + min (i , rand [ t +1] -v [k ]) [ imsd ])
116
- v [ imsd ] ) /
min (dv ,v[ imsu ] -v
; # print ( xi )
xw ← (v [k] - min (w ,v[ k] -Vmin ) ])
- v[ wmsd ] ) /
min (dv ,v [ wmsu ] -v [ wmsd
; # print ( xw )
117 118
Decisionmatrix ← cbind ( min (w ,v [k] -Vmin ) * ( P[,t ] *(1 - aw ) - bw ) + xw * C[, wmsu ] + (1 -xw ) *C[ , wmsd ] , C[ ,k] , ( P [,t ] *(1 + ai ) + bi
) + xi * C[, imsu ]
- min (i , rand [t +1] -v [k ]) * + (1 -xi )*C[ , imsd ] )
119 120
l ← apply ( Decisionmatrix ,1 , which.max ) - 2
121
l0 ← l ==0
122
l1 ← l ==1
123
lm1 ← l == -1
124 125
VM [ lm1 ,1 , k]
←
min (w ,v[ k] -Vmin ) * ( P [ lm1 ,t ] *(1 - aw ) - bw )
126
VM [ l0 ,1 , k ]
←
0
127
VM [ l1 ,1 , k ]
←
- min (i , rand [t +1] -v [k ]) * ( P[ l1 , t] *(1 + ai ) + bi
) 128
for (b in 1: Nsim ) {
129
VMtemp [b ,1 , k] ← ( xw *( lm1 [b ]) + xi *( l1 [b ]) ) * VM [b ,2 , k + l [b ]*( lm1 [ b ]) *( k-wmsu ) + l[ b ]*( l1 [b ]) *( imsu-k )]
130
+
((1 -xw ) *( lm1 [b ]) +(1 -xi ) *( l1 [b ]) ) * VM [b ,2 , k + l[ b
151
]*( lm1 [b ]) *( k-wmsd ) + l[b ]*( l1 [ b ]) *( imsd-k )] + ( l0 [b ]) * VM [b ,2 , k ]
131 132
}}
133
if ( rand [t +1]
134
wmsd ← k-nyceiling ( min (w , v[k ] -Vmin )/ dv )
135
wmsu ← wmsd + 1
136
xw ← (v [k] - min (w ,v[ k] -Vmin ) ])
;#
- v[ wmsd ] ) /
min (dv ,v [ wmsu ] -v [ wmsd
print ( xw )
137 138
Decisionmatrix ← cbind ( min (w ,v [k] -Vmin ) * ( P[,t ] *(1 - aw ) - bw ) + xw * C[, wmsu ] + (1 -xw ) *C[ , wmsd ] ,
(v [k] -rand [ t +1]) * ( P [,t] *(1 -
aw ) - bw ) + C[, which ( v == rand [t +1]) ] ) 139 140
l ← apply ( Decisionmatrix ,1 , which.max ) - 2
141
l0 ← l ==0
142
lm1 ← l == -1
143 144
VM [ lm1 ,1 , k]
←
min (w ,v[ k] -Vmin ) * ( P [ lm1 ,t ] *(1 - aw ) - bw )
145
VM [ l0 ,1 , k ]
←
146
for (b in 1: Nsim ) {
147
VMtemp [b ,1 , k] ← VM [b ,2 , k+ l[b ]*( lm1 [b ]) * min (w ,v[ k] -Vmin )/ dv - ( l0 [b
(v [k] -rand [ t +1]) * ( P[ l0 , t] *(1 - aw ) - bw )
]) *( v[ k] -rand [t +1]) / dv ] 148
}}}
149 150
# ######### Sidste k
151
if ( v[ Nm ]== rand [ t ]) {
152
k= Nm
153 154
wmsd ← k - 1 - nyceiling (( w- (v[ k] -v [ k-1 ]) )/ dv
152
)
155
wmsu ← wmsd + 1
156
xw ← (v [k] - min (w ,v[ k] -Vmin )
- v[ wmsd ] ) /
dv
157 158
Decisionmatrix ← cbind ( min (w ,v [k] -Vmin ) * ( P[,t ] *(1 - aw ) - bw ) + xw * C[, wmsu ] + (1 -xw ) *C[ , wmsd ] , (v [k] -rand [ t +1]) * ( P[ ,t] *(1 aw ) - bw ) + C[, which ( v == rand [t +1]) ])
159 160
l ← apply ( Decisionmatrix ,1 , which.max ) - 2
161
l0 ← l ==0
162
lm1 ← l == -1
163 164
VM [ lm1 ,1 , k]
←
min (w ,v[ k] -Vmin ) * ( P [ lm1 ,t ] *(1 - aw ) - bw )
165
VM [ l0 ,1 , k ]
←
166
for (b in 1: Nsim ) {
167
VMtemp [b ,1 , k] ← xw * VM [b ,2 , k+
(v [k] -rand [ t +1]) * ( P[ l0 , t] *(1 - aw ) - bw )
l[ b ]*( lm1 [b ]) *( k-wmsu ) - ( l0 [b ]) *( k-which (v ==
168
rand [ t +1]) )] + (1 -xw )* VM [b ,2 , k+ l[ b ]*( lm1 [b ]) *( k-wmsd ) - ( l0 [b ]) *( k-which (v ==
169
rand [ t +1]) )] 170
}}
171 172
# Updating step
173
VM [ ,2 , complete.cases ( VMtemp [1 ,1 ,]) ]← VMtemp [ ,1 , complete.cases ( VMtemp [1 ,1 ,]) ]* d+ VM [ ,1 , complete.cases ( VMtemp [1 ,1 ,]) ]
174
VM [ ,2 , is.na ( VMtemp [1 ,1 ,]) ] ← VM [ ,2 , is.na ( VMtemp [1 ,1 ,]) ]* d + w * P[ , t]
175
VM [ ,1 ,]← NA
176
VMtemp [ ,1 ,]← NA
177
Timemeasure [t]← proc.time () [ " elapsed "] -ptm [" elapsed "]
153
178
}
Listing 15: LMSC Hedge 2014-2015
Litteratur Nelson Areal, Artur Rodrigues, and Manuel R Armada. On improving the least squares monte carlo option valuation method.
Review of Derivatives Research, 11(1-2):119151,
2008. Fred Espen Benth, Steen Koekkebakker, and Fridthjof Ollmar. Extracting and applying smooth forward curves from average-based commodity contracts with seasonal variation.
The Journal of Derivatives, 15(1):5266, 2007. Alexander Boogert and Cyriel De Jong. Gas storage valuation using a monte carlo method.
The journal of derivatives, 15(3):8198, 2008. Alexander Boogert and Cyriel de Jong. Gas storage valuation using a multifactor price process.
The Journal of Energy Markets, 4(4):2952, 2011. The American Economic Review, pages 5072,
Michael J Brennan. The supply of storage. 1958.
Michael J Brennan and Eduardo S Schwartz. Evaluating natural resource investments.
Journal of business, pages 135157, 1985. John Breslin, Les Clewlow, Tobias Elbert, Calvin Kwok, and Chris Strickland. Gas storage: Overview & static valuation.
Lacima Group, 2008.
John Breslin, Les Clewlow, Tobias Elbert, Calvin Kwok, Chris Strickland, and Daniel van der Zee. Gas storage: Rolling instrinsic valuation.
Energy Risk, January, pages 6165, 2009.
Álvaro Cartea, James Cheeseman, and Sebastian Jaimungal. Natural gas storage modelling.
Handbook of Multi-Commodity Markets and Products: Structuring, Trading and Risk
Management, pages 877899, 2015. 154
Gonzalo Cortazar and Eduardo S Schwartz. Implementing a stochastic model for oil futures prices.
Energy Economics, 25(3):215238, 2003.
Rob J Hyndman.
forecast: Forecasting functions for time series and linear models.
Nicholas Kaldor. Speculation and economic stability.
The Review of Economic Studies, 7
(1):127, 1939. Nicholas Kaldor. A note on the theory of the forward market.
The Review of Economic
Studies, 7(3):196201, 1940. John Maynard Keynes. Treatise on money: Pure theory of money vol. i. 1930. Thomas Kremser and Margarethe Rammerstorfer. The convenience yield implied in the european natural gas markets-the impact of macrofactors and weather.
Available at SSRN
2146639, 2012. Francis A Longsta and Eduardo S Schwartz. Valuing american options by simulation: A simple least-squares approach.
Review of Financial studies, 14(1):113147, 2001.
Jan Müller, Guido Hirsch, and Alfred Müller. Modeling the price of natural gas with temperature and oil price as exogenous factors. In
Innovations in Quantitative Risk
Management, pages 109128. Springer, 2015. R Core Team.
R: A Language and Environment for Statistical Computing. R Foundation
for Statistical Computing, Vienna, Austria, 2014. URL http://www.R-project.org/. Paul A Samuelson. Proof that properly anticipated prices uctuate randomly.
Industrial
management review, 6(2):4149, 1965. Eduardo Schwartz and James E Smith. Short-term variations and long-term dynamics in commodity prices.
Management Science, 46(7):893911, 2000.
Eduardo S Schwartz. The stochastic behavior of commodity prices: Implications for valuation and hedging.
The Journal of Finance, 52(3):923973, 1997. 155
Lars Stentoft. Assessing the least squares monte-carlo approach to american option valuation.
Review of Derivatives research, 7(2):129168, 2004.
Sven-Olaf Stoll and Klaus Wiebauer. A spot price model for natural gas considering temperature as an exogenous factor and applications.
Journal of Energy Markets, 3(3):113128,
2010. Lester G Telser. Futures trading and the storage of cotton and wheat.
The Journal of
Political Economy, pages 233255, 1958. Thomas Volmer. Spot price models for natural gas. 2011. Holbrook Working. The theory of price of storage.
The American Economic Review, pages
12541262, 1949. Diethelm Wuertz, Yohan Chalabi with contribution from Michal Miklovic, Chris Boudt, Pierre Chausse, and others.
fGarch: Rmetrics - Autoregressive Conditional Heteroske-
dastic Modelling, 2013. URL http://CRAN.R-project.org/package=fGarch. R package version 3010.82.
156
