Niet-evenwichts moleculaire dynamica als een model voor financiële markten

Niet-evenwichts moleculaire dynamica als een model voor financiële markten Sander Van Goethem

Promotor: prof. dr. Jan Ryckebusch Masterproef ingediend tot het behalen van de academische graad van Master in de ingenieurswetenschappen: toegepaste natuurkunde

Vakgroep Fysica en Sterrenkunde Voorzitter: prof. dr. Dirk Ryckbosch Faculteit Ingenieurswetenschappen en Architectuur Academiejaar 2011-2012

Voorwoord Ik zou graag professor Jan Ryckebusch willen bedanken voor het onder de aandacht brengen van dit interessante onderwerp tijdens de lessen statistische fysica 2 en uiteraard ook voor de vele nuttige bedenkingen en tips tijdens het jaar. En niet vergeten het grote aantal spellingsen taalcorrecties. De auteur geeft de toelating deze masterproef voor consultatie beschikbaar te stellen en delen van de masterproef te kopiëren voor persoonlijk gebruik. Elk ander gebruik valt onder de beperkingen van het auteursrecht, in het bijzonder met betrekking tot de verplichting de bron uitdrukkelijk te vermelden bij het aanhalen van resultaten uit deze masterproef. The author gives permission to make this master dissertation available for consultation and to copy parts of this master dissertation for personal use. In the case of any other use, the limitations of the copyright have to be respected, in particular with regard to the obligation to state expressly the source when quoting results from this master dissertation. Ondertekend Sander Van Goethem, 4 juni 2012

i

Overzicht Sander Van Goethem Masterproef ingediend tot het behalen van de academische graad van Master in de ingenieurswetenschappen: toegepaste natuurkunde Titel: Niet-evenwichts moleculaire dynamica als een model voor financiële markten Promotor: prof. dr. Jan Ryckebusch Faculteit Ingenieurswetenschappen en Architectuur Vakgroep Fysica en Sterrenkunde Voorzitter: prof. dr. Dirk Ryckbosch Faculteit Ingenieurswetenschappen en Architectuur Academiejaar 2011-2012

Samenvatting We proberen in deze masterproef de dynamica van de financiële markten na te bootsen met behulp van Niet Evenwichts Moleculaire Dynamica (NEMD). Met behulp van de returns van de DAX index in de periode 1991-2011 wordt de dynamica onderzocht. Er worden vette staarten vastgesteld in de probabiliteitsdistributie functie (PDF) van deze returns. Een link wordt gelegd tussen het economische begrip volatiliteit en de fysische informatiëentropie. Deze eigenschappen van de Duitse DAX index worden vergeleken met een zeer eenvoudig economisch model gebaseerd op de geometrische brownse beweging (GBB). We stellen vast dat dit GBB model faalt in het beschrijven van de vette staarten die een machtswet vertonen en dat daarom een nieuw soort model nodig is. Dit model gaan we construeren in het raamwerk van het goed gekende klassieke moleculaire dynamica (MD) model. Hierna volgt een korte beschrijving van dit model en een theoretische uitweiding over de entropie. Er wordt aangetoond dat de returns van de marktindexen kunnen gerelateerd worden aan de stapgroottes van de deeltjes in de MD simulatie. In hoofdstuk 5 wordt dan geprobeerd de MD simulatie op een bepaalde manier uit evenwicht te brengen zodat de stapgroottes van de deeltjes een ii

iii

PDF verkrijgen met staart die aflopen zoals machtswetten. Met veel ’trail and error’ wordt er uiteindelijk een PDF verkregen met een machtswetstaart. Voorzichtigheid is echter geboden aangezien het resultaat een toevalstreffer kan zijn. Verder onderzoek is nodig.

Trefwoorden econofysica, niet-evenwichts moleculaire dynamica, machtswetstaarten

Non-equilibrium molecular dynamics as a model for the financial markets Sander Van Goethem Supervisor(s): prof. dr Jan Ryckebusch Abstract— The dynamics of the financial markets are examined by analyzing the daily closing time data of the German DAX index in the period 1991-2011. We then try to replicate the characteristic powerlawtails observed in the probabilitydistributionfunction (PDF) of the returns of the DAX index with a molecular dynamics (MD) simulation which is being driven out of equilibrium by tinkering with the temperature of the system. Keywords—Econophysics, non equilibrium molecular dynamics (NEMD), powerlawtails

Fig. 1. The positive and negative cumulative distribution of the returns of the DAX index. To both functions a powerlaw (y = Cxα ) is fitted.

I. I NTRODUCTION UE to the large impact and importance of certain market phenomena, such as crashes and periods of high volatility, there has been a lot of effort to try and model the complex dynamics of the markets in a mathematical way. In recent times physicists have also excepted this challenge ([1],[2]), naming a new branch of physics: econophysics. In this article we try to establish a connections between the financial system and a classical fluid simulated with NEMD.

D

II. A NALYSIS OF THE DAX INDEX First we analyze the properties of the DAX index in the period 1991-2011 to figure out the dynamics of the index. A. Returns The logarithmic return S(t) is defined as Y (t + ∆t) . (1) Y (t) with Y(t) the value of the index at time t. We choose ∆t to be equal to 1 day. Figure 1 shows that the CDF has powerlawtails which implies that also the PDF has tails that show powerlaw behavior [3]. An exponentially decaying autocorrelation function indicate that the returns of the DAX index are short-ranged correlated. S(t) = ln

B. Volatility/Informationentropy Volatility is a measure that quantifies risk and is an important quantity for investors. The volatility over a certain period T = n∆t is calculated from the returns S(t) as follows Pτ (S(i) − µ)2 σT (τ ) = i=τ −T . (2) T −1 We can also define a time-dependent informationentropy, calculated out of the PDF p of the returns of the DAX index. The summation runs over the bins b. H(τ ) =

X b

−pb [S(τ − T ) → S(τ )] ln pb [S(τ − T ) → S(τ )] . (3)

For a properly chosen value of T the volatility and the entropy show a similar dynamic behavior. In figure 2 both are shown for T = 200. This unique similarity is highly dependent on the chose of T . For T = 30, for example, the similarity is gone.

Fig. 2. The volatility (red) and informationentropy (green) with a period of T = 200 days of the DAX index between 1991 and 2011.

III. M OLECULAR DYNAMICS MD is a well known technique to simulate classical molecular systems. The Newtonian equations of motion are integrated numerically at regular time intervals for many interacting particles with the Velocity Verlet Algorithm. We execute the simulation with a system consisting of 256 argon atoms and simulate a period of 100, 000 timesteps. A. Softcore Potential The forces that act on the particles are derived from an empirical potential. In equilibrium MD it is usual to use the LennardJones potential which incorporates the Pauli repulsion at short ranges and the long-range van der Waals force. When the system is brought out of equilibrium it is possible for particles to approach very close to each other. Because of the high values of the LJ potential in these cases the simulation can become unstable. In order to remedy this, the LJ potential is replaced by the softcore (SC) potential [4].

its original value. To assure that the temperature remains within bounds, the temperature change is not induced if the random number is too high. B. Result The result of these temperature manipulations are shown in figure 3. We plotted the positive tails of the CDF’s of the returns of the DAX index S(t) and of the position displacement ∆x(t). We observe that the proposed temperature changes drive the system of argon atoms in such a way that a powerlawtail, that is similar to the DAX index, is generated. Fig. 3. The distribution function of the positive tails of both the returns of the DAX index (green) and the proposed NEMD simulation (blue).

(r − RA )2 H − UA exp − . 2 1 + exp(∆(r − RR )) 2δA (4) The parameters in the SC potential are optimized to match the long-range behavior of the LJ potential. H is called the hardness parameter and determines the height of the potential at the core. USC (r) =

B. Link with markets At every timestep in the simulation the average displacement of all the particles ∆x is saved. The return of the DAX index S(t) represents a change in value of the DAX index, the average displacement ∆x(t) respresents a change in position of the particles in the MD simulation. When properly normalized the two can be related to each other. ∆x(t) − h∆x(t)i S(t) − hS(t)i ←→ . σS(t) σ∆x(t)

(5)

IV. T EMPERATURE DRIVEN NEMD The chosen method of driving the system out of equilibrium consists of abruptly changing the temperature of the system. The amplitude of the temperature change is determined by random generated numbers. A. Inducing temperatureshocks Every 1, 000 steps a random number rnd is generated. 70% of the time it is extracted from a gaussian distribution with µ = 0 and σ = 1, and 30% of the time it is extracted from a Lévy alpha-stable (LAS) distribution, defined by its characteristic function given below, with parameters a = 0.01, β = 0, µ = 0.5 and m = 0. h i ˆ a,β,m,µ (z) = exp −a|z|µ 1 + iβsgn(t) tan( πµ ) + imz . L 2 (6) When the random number rnd is positive, the temperature is multiplied by (1 + rnd), when its negative it is divided by (1 − rnd). After another 200 steps the temperature is restored to

V. C ONCLUSION The proposed NEMD simulation generates powerlawtails in the PDF of the position displacement, a robust feature of the PDF of the returns of marketindices. However, a lot of caution has to be exercised. Because of the randomness of the temperaturechange and the limited simulation time, it is possible that the resulting outcome was a ’lucky hit’. More and longer simulations are necessary before we can call this simulation a succes. R EFERENCES [1] J. Voit, Statistical Mechanics of Financial Markets, Springer, 3th edition, 2005. [2] H. Eugene Stanley Rosario N. Mantegna, Introduction to econophysics, Cambridge University Press, 2005. [3] M. E. J. Newman, “Pareto laws, pareto distributions and zipf’s law,” Contemporary Physics, vol. 46, pp. 323–351, 2005. [4] J. Ryckebusch S. Standaerd and L. De Cruz, “Creating conditions of anomalous self-diffusion in a liquid with molecular dynamics.,” Journal of Statistical Mechanics, p. P04004, 2010.

Inhoudsopgave Voorwoord

i

Overzicht

ii

Inhoudsopgave

vi

Lijst van Afkortingen en Symbolen

ix

1 Inleiding

1

2 Statistiek van financi¨ ele markten

2

2.1

2.2

Begrippen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2

2.1.1

Efficient-market hypothesis . . . . . . . . . . . . . . . . . . . . . . . .

2

2.1.2

Returns . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3

2.1.3

Volatiliteit

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4

2.1.4

Autocorrelatiefunctie . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6

2.1.5

Stochastische processen . . . . . . . . . . . . . . . . . . . . . . . . . .

6

2.1.5.1

Martingale proces . . . . . . . . . . . . . . . . . . . . . . . .

6

2.1.5.2

Markov proces . . . . . . . . . . . . . . . . . . . . . . . . . .

6

2.1.5.3

Wiener proces . . . . . . . . . . . . . . . . . . . . . . . . . .

7

2.1.6

Centrale Limiet Theorema (CLT) . . . . . . . . . . . . . . . . . . . . .

7

2.1.7

Lévy alpha-stable distributie (LAS) . . . . . . . . . . . . . . . . . . .

8

2.1.8

Informatieëntropie . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8

Geometrische Brownse Beweging . . . . . . . . . . . . . . . . . . . . . . . . .

11

2.2.1

Tijdschaal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

13

2.2.2

Returns . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

15

2.2.3

Volatiliteit

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

17

2.2.4

Informatieëntropie . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

19

2.2.5

Autocorrelatie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

22

vi

vii

Inhoudsopgave

2.2.6

Conclusie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3 Moleculaire Dynamica 3.1

3.2

23 25

Oplossingsmethode . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

26

3.1.1

Verlet algoritme . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

26

3.1.2

Interactiepotentiaal

. . . . . . . . . . . . . . . . . . . . . . . . . . . .

26

3.1.3

Systeemheden . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

27

Programma . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

29

3.2.1

Initialisatie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

29

3.2.2

Thermalisatie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

29

3.2.3

Productiefase . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

30

3.2.3.1

Temperatuur . . . . . . . . . . . . . . . . . . . . . . . . . . .

30

3.2.3.2

Energie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

31

3.2.3.3

∆x . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

32

3.2.3.4

Entropie . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

34

3.2.3.5

Volatiliteit . . . . . . . . . . . . . . . . . . . . . . . . . . . .

36

3.2.3.6

Snelheidsautocorrelatiefunctie . . . . . . . . . . . . . . . . .

36

4 Entropie

38

4.1

Entropie in de Statistische Mechanica . . . . . . . . . . . . . . . . . . . . . .

39

4.2

Entropie van veeldeeltjessysteem . . . . . . . . . . . . . . . . . . . . . . . . .

41

4.2.1

Translatieëntropie . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

41

4.2.2

Interactieëntropie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

43

Entropie in de Simulatie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

45

4.3

5 Simuleren van Markten met NEMD

46

5.1

Temperatuurschok . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

46

5.2

Temperatuursgradiënt . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

49

5.3

Marktsimulatie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

52

5.3.1

’Gaussische’ schokken . . . . . . . . . . . . . . . . . . . . . . . . . . .

52

5.3.2

’Levy’ schokken . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

57

5.4

’Levy-Gauss’ schokken . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

60

5.5

Opmerkingen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

62

6 Conclusie

65

Bibliografie

66

Inhoudsopgave

viii

Lijst van figuren

68

Lijst van tabellen

73

Lijst van Afkortingen en Symbolen BEL20

Belgisch beursindex

C(τ )

autocorrelatiefunctie

CDF

Cumulatieve Distributiefunctie

CLT

Centrale Limiet Theorema

DAX

Deutscher Aktienindex, de Duitse beursindex

E

energie

Ek

kinetische energie

Ep

potentiële energie

EHM

efficiënt markt hypothesis

GBB

Geometrische Brownse Beweging

h

constante van Planck

H

informatieëntropie

H

hardheidsparameter in de SC potentiaal

kb

constante van Boltzmann

LAS

Lèvy alpha-stable distributie

LJ

Lennard-Jones (potentiaal)

m

massa

MD

Moleculaire Dynamica

N

aantal deeltjes

NEMD

niet-evenwichts moleculaire dynamica

NIKKEI

Japans beursindex

p

momentum

PDF

Probabiliteitsdistributiefunctie

P(x)

Probabiliteitsdistributiefunctie

S

Entropie

S

Returns

SDV

Statistische Differentiaalvergelijking

SP500

Standaard & Poor’s 500, beursindex van de 500 grootste verhandelde aandelen in de VS ix

x

Inhoudsopgave

SC

Softcore (potentiaal)

t

tijd

T

Temperatuur

U

interactiepotentiaal

Y

prijs van een asset

Z

’Zustandsumme’, partitiefunctie

β

1 kb T

energieschaal in de LJ potentiaal

µ

gemiddelde

λDB

de Broglie golflengte

ρ

dichtheid

σ

standaardafwijking

σ

lengteschaal in de LJ potentiaal

σt

volatiliteit (per t dagen)

Hoofdstuk 1

Inleiding Ingenieurs en natuurkundigen gebruiken in grote mate statistische methoden om complexe processen die niet kunnen beschreven worden door analytische methoden beter te begrijpen. Sinds geruime tijd wordt deze statistische en kritische blik ook gericht op een van de meest complexe systemen gecreëerd door de mens: de economie en in het bijzonder de grillen van de beurs. Het onderzoeken van de markten met behulp van methoden uit de statistische fysica wordt sinds 1995 dankzij E. Stanley ondergebracht onder de term ’econofysica’. De statistische beschrijving van de markt is de gemakkelijkste wiskundige manier om de interactie tussen miljoenen beleggers, speculanten en ’hedgers’ te beschrijven. In deze masterproef wordt gepoogd de complexe dynamica van de financiële markten na te bootsen met behulp van het gekende klassieke moleculaire dynamica model. Eerst wordt er marktdata onderzocht om de dynamica van de markt bloot te leggen. Door de vrij verkrijgbare daily closing time data [1] over de DAX index in de periode 1991-2011 te analyseren met statistische methoden, verkrijgen we bepaalde grootheden, met name de volatiliteit en de informatieëntropie van de returns, die we gaan proberen linken aan moleculaire grootheden uit de MD simulatie. Hierna volgt een korte beschrijving van de gebruikte MD simulatie en een theoretische behandeling over de entropie. Vervolgens worden de vergaarde inzichten gebruikt om het fysische systeem in de MD simulatie op een bepaalde manier uit evenwicht te brengen. Dit gebeurt door het systeem te onderwerpen aan plotse temperatuursveranderingen. Zo verkrijgen we anomale diffusie die we kunnen linken aan de extreme events vanop de financiële markten. Uiteindelijk worden uit de bevindingen conclusies getrokken over bekomen resultaten en de juistheid en relevantie ervan.

1

Hoofdstuk 2

Statistiek van financi¨ ele markten De wiskundige analyse van de financiële markten kent zijn oorsprong in de beroemde paper van Louis Bachelier ’La Theorie de la spéculation’ [2]. Vanaf dan is er zeer veel tijd en moeite ge¨ınvesteerd in het bepalen van de onderliggende waarde van aandelen, futures en opties via mathematische modellen. Een voorbeeld is het veelgebruikte Black-Scholes model voor de prijs van opties, beschreven door Fischer Black en Myron Scholes in de vroege jaren zeventig [3]. In dit hoofdstuk zullen we een eenvoudig model bespreken waarvan aangenomen wordt dat het de dynamica van de markten redelijk beschrijft. Na het introduceren van een paar begrippen, zullen we de Geometrische Brownse Beweging bespreken en de link leggen met de financiële markten. Vervolgens worden de statistische variabelen die uit dit model berekend kunnen worden, vergeleken met die van de Duitse DAX index in de periode 1991-2011.

2.1

Begrippen

Voor een goed begrip zullen eerst een paar concepten en definities worden gegeven die we later nog zullen gebruiken.

2.1.1

Efficient-market hypothesis

De efficient-market hypothesis (EMH) stelt dat de financiële markt op elk moment alle relevante informatie ter beschikking heeft en deze instantaan en rationeel gebruikt. Hierdoor zal de prijs van een asset op elk moment de ideale prijs zijn (de echte waarde van het asset, de waarde die de behoeften van de investeerders reflecteert) die alle beschikbare informatie uit het verleden en het heden in rekening brengt. De evolutie van de prijs van het asset zal dus enkel afhangen van de nieuwe informatie die de markt binnensijpelt. De 2

Hoofdstuk 2. Statistiek van financiële markten

3

ef f icientmarkethypothesis is uiteraard maar een benadering en er is discussie over hoe goed deze wel is en op welke punten een markt zich anders gedraagt.

2.1.2

Returns

Een belangrijke variabele voor het analyseren van de markten is de return. Voor een belegger is het belangrijk dat men een idee heeft van de prestaties van zijn investering. De relatieve stijging of daling van de marktprijs van een asset kan worden samengevat in de returns. Om de returns te berekenen van de prijsevolutie van een asset zijn er verschillende methoden. Zij Y (t) de prijs van een bepaald asset op tijd t, dan kunnen we de return van de prijs Y (t) definiëren op de volgende manieren. • een eerste eenvoudige manier om de return te definiëren is het verschil te nemen van de prijs op opeenvolgende tijdstippen

Sversch (t) = Y (t + ∆t) − Y (t) , • we kunnen ook een procentuele return definiëren als Sproc (t) =

Y (t + ∆t) − Y (t) , Y (t)

• als laatste is er ook nog de logaritmische return Y (t + ∆t) Slog (t) = ln . Y (t) Merk op dat voor zeer kleine ∆t de procentuele en de logaritmische return samenvallen. In het volgende zullen we enkel gebruik maken van de logaritmische return. In figuur 2.1 is een voorbeeld gegeven van de genormaliseerde logaritmische returns van de BEL20 index in de periode van juli 2005 tot begin 2011. De return is berekend voor een ∆t van 1 dag en met de waarde van de index op de sluitingstijd van de beurs. Op deze manier zullen de returns in deze masterproef altijd berekend worden. Vooral interessant voor de beleggers en diegene die de markt correct willen simuleren, zijn de grote uitschieters die hoge winsten en verliezen voorstellen. Om de relevante data van de verschillende indexen en de simulaties te kunnen vergelijken, zullen telkens de returns genormaliseerd worden. Dit kan bewerkstelligd worden door de gemiddelde en de variantie te berekenen over een bepaalde periode T . µ=

P

S(t) , N

T

(2.1)

4


Figuur 2.1: De logaritmische returns (∆t = 1 dag) van de BEL20 index in de periode van juli 2005 tot begin 2011. Om de returns te berekenen is de slotkoers van elke beursdag gebruikt [1]. De returns vertonen op bepaalde momenten uitgesproken pieken die overeenkomen met grote winsten en verliezen (zie bv de crash in september 2008). Deze pieken, en het feit dat ze geclusterd zijn, worden ook voor bijna elke andere index teruggevonden.

P ( T S(t) − µ)2 σ = . N −1 De genormaliseerde returns worden dan gegeven door 2

Sn (t) =

2.1.3

S(t) − µ . σ

(2.2)

(2.3)

Volatiliteit

De standaard manier om het risico van een bepaalde investering te bepalen is aan de hand van de volatiliteit. Dit is gedefinieerd als de standaardafwijking van de returns over een bepaald tijdsinterval T .

5


σT (τ ) =

Pτ

i=τ −T (S(i)

T −1

− µ)2

.

(2.4)

Om een statistisch significante volatiliteit uit te komen moeten we het tijdsinterval groot genoeg kiezen. Hoe groter de volatiliteit, hoe groter de kans op grote verliezen en winsten. Op die manier kan de volatiliteit dan ook gezien worden als een maat voor de nervositeit op de financiële markten, een nervositeit die ervoor zorgt dat bv. de EMH niet meer opgaat en de markt zich alles behalve rationeel gedraagt. In figuur 2.2 is een voorbeeld gegeven van de volatiliteit van de BEL20 index voor T = 30 beursdagen, berekend met returns van daily closing time data. De grote piek in de tweede helft van 2008 is een gevolg van de grote financiële crisis die op dat moment de kop op stak. Figuur 2.2: Volatiliteit van de BEL20 index in de periode van juli 2005 tot begin 2011. De volatiliteit kan gezien worden als de nervositeit of de ’temperatuur’ van de markten. Dit is duidelijk te zien in het tweede deel van het jaar 2008 door de financiële crisis in de tweede helft van 2008.


2.1.4

6

Autocorrelatiefunctie

Om een beeld te krijgen van de tijdscorrelaties in een systeem kan men gebruik maken van de autocorrelatiefunctie. Deze kan op verschillende manieren gedefinieerd worden, maar in dit werk zullen we enkel gebruik maken van de volgende definitie. C(τ ) =

hX(t + τ )X(τ )i − hX(t + τ )i ∗ hX(t)i . σ2

(2.5)

C(τ ) geeft weer hoe de waarde van een variabele X(t) de waarde van de variabele X(t + τ ) be¨ınvloedt. De tijdscorrelaties van de statistische variabelen kunnen meer inzicht bieden in het proces dat zich op de financiële markten afspeelt.

2.1.5

Stochastische processen

Een stochastisch proces is een model voor een proces met als uitkomst stochastische variabelen. Met andere woorden processen waarvan de uitkomst van het toeval afhangt. De schommelingen op de beurzen worden als een stochastisch proces beschouwd. Een paar veelgebruikte stochastische processen zijn de volgende. 2.1.5.1

Martingale proces

Een discrete tijd martingaal proces is een discrete tijd stochastisch proces met X een random variabele dat op elk tijdstip t aan de volgende voorwaarde voldoet E(Xt+1 |X1 , X2 , .., Xt ) = Xt , of equivalent E(Xt+1 − Xt |X1 , X2 , .., Xt ) = 0 . E(Xt+1 − Xt |X1 , X2 .., Xt ) wordt de conditionele verwachtingswaarde genoemd en wordt berekend met behulp van de conditionele probabiliteit p(Xt+1 − Xt |X1 , X2 , .., Xt ) die de kans is dat (Xt+1 − Xt ) zich voordoet als (X1 ,X2 , .. ,Xt ) gekend zijn. 2.1.5.2

Markov proces

Een Markov proces is een proces waarbij de waarde van de random variable op tijdstip t + 1 Xt+1 enkel afhangt van de waarde van de random variable op tijdstip t Xt . p(Xt+1 |Xt , Xt−1 , ..., X1 ) = p(Xt+1 |Xt ) .

(2.6)

7


Dit wil zeggen dat het proces geen geheugen heeft: alle informatie over het verleden (Xt−1 ,Xt−2 ,...,X1 ) is irrelevant voor de toekomst. We zullen dit vaak terug vinden in de vooropgestelde modellen voor de dynamica van markten. 2.1.5.3

Wiener proces

Een speciaal Markov proces is het Wiener proces. Dit proces, beroemd geworden door Einstein’s beschrijving van de Brownse beweging in 1905 [4], was eigenlijk al voorgesteld door Bachelier in 1900 [2]. Een speciale eigenschap van het Weiner proces is dat de verandering van de stochastische variabele (Xt+1 − Xt ) gegeven wordt door √ (Xt+1 − Xt ) = ∆t .

(2.7)

Voor een voldoende kleine waarde van ∆t = (t + 1) − t en met getrokken uit een Gaussische

distributie met µ = 0 en σ = 1.

2.1.6

Centrale Limiet Theorema (CLT)

Stel dat Sn de som is van n onafhankelijke randomvariabelen xi , met E(xi ) = 0 en E(x2i ) = s2i (eindige variantie) Sn =

n X

xi ,

(2.8)

i=1

met

σn = E(Sn2 ) =

n X

s2i .

(2.9)

i=1

Dan

1

zegt het centrale limiet theorema dat Sn Sen = , σn

Gaussisch verdeeld is met E(Sen ) = 0 en E(Sen2 ) = 1

P Sen

1 = √ exp 2π

(2.10)

Sen2 2

!

.

(2.11)

Als aangenomen wordt dat de returns van een index randomvariabelen zijn met een eindige variantie dan verwachten we dat de returns van een financiële index normaal verdeeld zijn. 1

Er moet ook nog voldaan worden aan de Lindeberg conditie [5]


2.1.7

8

L´ evy alpha-stable distributie (LAS)

De Lévy alpha-stable distributies worden gedefinieerd door de volgende karakteristieke vergelijking: h i ˆ a,β,m,µ (z) = exp −a|z|µ 1 + iβsgn(t) tan( πµ ) + imz . L 2

(2.12)

De parameter a wordt de schaalparameter genoemd en bepaalt de breedte van de distributie. β is een maat voor de gepiektheid en de symmetrie van de distributie. De locatie van het centrum van de distributie wordt bepaald door m. µ is de stabiliteit van de distributie en bepaalt (als µ < 2) het gedrag van de distributie op oneindig. Een speciale eigenschap van deze distributies is dat ze stabiel zijn. Dit wil zeggen dat de som van random variabelen die verdeeld zijn volgens een LAS distributie, ook een LAS distributie vormt. Als Sn de som is van n onafhankelijke randomvariabelen xi , met E(xi ) = 0 en E(x2i ) = ∞ dan zal (voor µ < 2) Sn verdeeld zijn volgens een LAS distributie [6].

Een andere interessante eigenschap is dat de staarten van de LAS distributies een machtswet volgen. Veronderstel dat β = 0 (distributie symmetrisch) en m = 0 (distributie gecentreerd rond 0), dan is voor µ < 2 µ ˆ µ (x) ≈ µA , voor x → ∞ . L |x|1+µ

2.1.8

(2.13)

Informatie¨ entropie

Informatie-entropie is een begrip ingevoerd door Claude Shannon [7]. Ge¨ınspireerd door het fysische begrip entropie (waarover later meer), is informatieëntropie een maat voor de informatie gestockeerd in een reeks onafhankelijke randomvariabelen xi . Met p(xi ) de probaliteitsdistributiefunctie (PDF) van de randomvariabele xi wordt de informatie-entropie H in het geval van een continue variabele gegeven door: Hcont = −

Z

xmax

dx p(x) ln[α p(x)] ,

(2.14)

xmin

De parameter α moet ingevoerd worden om er voor te zorgen dat de dimensie juist is. Er geldt immers dat Z

xmax

dx p(x) = 1 ,

(2.15)

xmin

waaruit volgt dat [p(x)] =

1 [x] .

Daar de entropie een dimensieloze grootheid is moet −ln[α p(x)]

ook dimensieloos zijn. Om dit dimensieloos te maken is er dus een α nodig met dimensie [x].

9


Er wordt dan ook gezegd dat α de karakteristieke lengteschaal in het probleem is [8]. Deze α is vrij te kiezen en bepaalt mee de absolute grootte van de entropie. De financiële datareeksen die hier onderzocht worden zijn echter bij definitie discreet. In het discrete geval kan vergelijking 2.14 benaderd worden, met N het totaal aantal bins en ∆ de bingrootte, door

Hdisc

∼ =− ∼ =−

N X

i=0 N X i=0

∆p(

xi−1 + xi xi−1 + xi ) ln p( ), 2 2

∆p(

X xi−1 + xi xi−1 + xi xi−1 + xi ) ln ∆p( ) + ln ∆ ). ∆p( 2 2 2

(2.16) N

(2.17)

i=0

Daar voor N >>> 1 N X

∆p(

i=0

vinden we

Hdisc − ln ∆ = −

N X

xi−1 + xi )=1, 2

∆p(

i=0

xi−1 + xi xi−1 + xi ) ln ∆p( ). 2 2

(2.18)

(2.19)

De bingrootte zal dus de waarde van de entropie be¨ınvloeden. Wederom wordt de entropie be¨ınvloed door een lengteparameter, in dit geval de bingrootte. Dit is te zien in figuur 2.3 waar we de entropie van de BEL20 index bereken aan de hand van verschillende bingroottes. De bingrootte ∆ is omgekeerd evenredig met het aantal bins. Wat in feite gedaan is, is dat we telkens de PDF p van de returns berekend hebben met een afnemende bingrootte. Uit deze PDF’s is dan telkens de informatieëntropie berekend volgens vergelijking 2.19. Om consistent te kunnen vergelijken moeten we dus alle datasets normaliseren en de bingrootte constant houden, wat in het continue geval neerkomt op het vastleggen van α. De informatieëntropie van de financiële markt zal een bepaalde waarde hebben afhankelijk van de distributie van de returns. Deze waarde op zich zegt ons niet veel. Daarom wordt de bekomen waarde telkens vergeleken met de waarde van de entropie voor een Gaussische distributie, uiteraard berekend met dezelfde bingrootte. Er kan bewezen worden [8] dat als Z

+∞

−∞

en

dxp(x) = 1 ,

(2.20)


10

Figuur 2.3: In deze figuur is de informatieëntropie (Hdisc − ln ∆) weergegeven in functie van het

aantal bins. De x-as is in logaritmische schaal. De entropie stijgt in functie van het aantal bins (∆ → 0)

Z

+∞

dxx2 p(x) = σ 2 ,

(2.21)

−∞

de entropie een maximale waarde heeft als x2 exp − 2 p(x) = √ 2σ 2πσ 2 1

(2.22)

Dit wil dus zeggen dat als de returns verdeeld zijn volgens een gaussische distributie de entropie maximaal is.


2.2

11

Geometrische Brownse Beweging

Het model dat we gaan gebruiken om te vergelijken met de reële markt, is het zogenaamde Îto proces (beschreven in bv. [8]) dat ook wel geometrische Brownse beweging (GBB) wordt genoemd. Dit model dat zijn wortels heeft in de paper van Bachelier en populair is geworden dankzij Einstein’s theorie over de geometrische Brownse beweging, wordt al bijna een eeuw gebruikt (ook vandaag de dag nog [9], [8]) en wordt gezien als een redelijk laagste orde model voor de financiële markten. De stochastische differentiaal vergelijking (SDV) is ∆Y (t) = µY (t)∆t + σY (t)∆X(t) .

(2.23)

De eerste term in deze SDV stelt dat de grootte van de verandering ∆Y van de assetprijs Y evenredig is met de assetprijs zelf. De term µ∆t wordt de drift term genoemd en wordt ingevoerd om de stijgende (of dalende) trend die de markten vertonen in het model te incorpo√ reren. De tweede term bevat het stochastisch proces (het Wiener proces 2.7) ∆X(t) = ∆t met getrokken uit een gaussische distributie met gemiddelde 0 en standaardafwijking 1. De waarde van σ bepaalt dan hoe sterk de prijs fluctueert. In figuur 2.4 tonen we een simulatie van het Îto proces met verschillende parameters. Het GGB model stelt een ’ideale’ financiële markt voor: constante groei (door de drift term µ) en redelijk stabiel (door de lage σ). Om te vergelijken met andere indexen kunnen we echter niet zomaar een willekeurige waarde voor µ en σ kiezen. De parameters µ en σ worden bepaald door het gemiddelde en de standaardafwijking te berekenen van de returns van de index. In figuur 2.5 wordt de prijsevoluties en de bij behorende GBB simulatie getoond van een index met positieve drift, de SP500. Een index met negatieve drift, de NIKKEI wordt afgebeeld in figuur 2.6. Een GBB simulatie is een random proces. Het is dus ook normaal dat de eindwaarde in de simulaties niet de eindwaarden van de index is. In het geval van de GBB simulatie van de NIKKEI is dit nu in deze figuur wel het geval maar in andere simulaties kan de waarde veel hoger of lager liggen zoals bv. bij de simulatie van de SP500 index het geval is. We kunnen de µ van de GBB simulatie ook vertalen in een jaarlijkse procentuele groei. In deze twee voorbeelden komen µ = 0, 0001565 en µ = −0, 000076213 overeen met een jaarlijkse groei in de GBB simulatie van 3, 44% en −4, 98%. De groei van de indexen zelf is resp. 6, 92% en −4, 97% per jaar in de periode 1991-2011. De procentuele groei in een GBB simulatie kan zeer sterk verschillen van simulatie tot simulatie en zegt niets over de correctheid waarmee een GBB simulatie een bepaalde index kan simuleren.


12

Figuur 2.4: De figuur toont een GBB simulatie met verschillende parameters (∆t = 0, 01 en 10000 tijdstappen). De parameter µ bepaalt hoe snel de prijs stijgt, σ bepaalt daarintegen hoe sterk de prijs fluctueert rond de exponentiële groei.

In de praktijk zijn zowel drift term µ als de σ niet constant maar tijdsafhankelijk wat een serieuze beperking oplegt aan de GBB simulatie. GGB kan veralgemeend worden in de zogenaamde ARCH en GARCH [10] processen die rekening houden met de tijdsafhankelijkheid van de twee parameters µ en σ. Deze uitbreidingen zullen hier niet behandeld worden. De parameters µ en σ die we in het GBB model moet stoppen kunnen we op twee verschillende manieren bepalen. Enerzijds kunnen we een Gaussische distributie fitten aan de PDF van de returns van de marktdata, anderzijds kunnen we gewoon het gemiddelde (µ) en de standaardafwijking (σ) van de returns berekenen. In figuur 2.7 worden de resultaten van deze twee verschillende methoden getoond. De rode PDF is de PDF van de returns van een index (in dit geval de DAX). De groene PDF is de Gaussische fit aan de PDF van de returns van de DAX index en de blauwe is een Gaussische functie met het gemiddelde gelijk aan het gemid-


13

Figuur 2.5: Een GBB simulatie met parameters bepaalt uit het gemiddelde en de standaardafwijking van de te simuleren SP500 index in de periode 1991-2012(µ = 0, 0001565, σ = 0, 0052).

delde van de returns van de DAX en de standaardafwijking gelijk aan de standaardafwijking van de DAX. Op figuur 2.7 is te zien dat de blauwe PDF de rode het ’best’ benadert. We zullen dus de parameters voor de GBB telkens bepalen door het gemiddelde en de standaardafwijking te berekenen uit de returns van de indexen die we willen simuleren. Merk ook nog op dat in vergelijking met een Gaussische distributie kleine returns meer waarschijnlijk zijn, middelmatige returns minder waarschijnlijk en grote returns veel waarschijnlijker.

2.2.1

Tijdschaal

Voor we de dynamica van de GBB simulatie beginnen te vergelijken met de dynamica van de financiële markt, is het belangrijk een paar opmerkingen te maken over de tijdschalen. Vandaag de dag wordt elke verrichting op de beurs opgevolgd en verandert de prijs van het asset instantaan. Deze hoog frequente data bevat veel informatie en wordt vaak gebruikt om statistische analyse op toe te passen. Het grote probleem met deze data is dat er een bepaalde Trading Time is. De beurzen zijn niet 24 op 24 open dus als we een periode van een paar dagen willen analyseren en we dus continue data nodig hebben, zullen we de data van de verschillende Trading dagen aan elkaar moeten plakken. De verandering van de prijs


14

Figuur 2.6: Een GBB simulatie met parameters bepaalt uit het gemiddelde en de standaardafwijking van de te simuleren NIKKEI index in de periode 1991-2012 (µ = −0, 000076213, σ =

0, 0067).

bij de sluiting en bij de opening de volgende dag wordt dus behandeld als een verandering op zeer korte tijd. Het is duidelijk dat gedurende de tijd dat de beurs gesloten is er toch nog informatie de markt binnensluipt uit bv. andere delen van de wereld en de verandering van de prijs bij de sluiting en bij de opening de volgende dag niet instantaan gebeurt zonder extra informatie. We kunnen langere perioden ook op een andere manier analyseren door bv. enkel de prijs bij sluiting te beschouwen. De tijdstap kunnen hier wel als constant beschouwen (eigenlijk niet: tijdens bv het weekend zijn de beurzen ook niet open) maar we zijn zeer veel informatie verloren. Een ander voordeel van het gebruik van de slotkoersen is het feit dat deze voor de belangrijkste indexen gemakkelijk vrij te verkrijgen zijn op het internet [1]. In deze thesis zullen we altijd de dagelijkse sluitingsprijzen analyseren. Het dynamisch gedrag van de financiële markt hangt ook af van de tijdschaal. Tegenwoordig wordt er bv. zeer veel aan zogenaamd ’algorithmic trading’ gedaan: het aangaan van transacties met behulp van computeralgoritmes die veel sneller kunnen handelen en beslissen dan mensen dat kunnen. Dit levert geheel ander gedrag op de subseconde tijdsschaal [11]. Het is duidelijk dat als we niets van de informatie op microseconde schaal behouden, deze dynamica


15

Figuur 2.7: Het bepalen van de parameters voor de GBB simulatie kan gebeuren door: 1. Een Gaussische functie te fitten (groen) aan de PDF van de returns (rood) of 2. Door met het gemiddelde en de standaarafwijking van de returns te berekenen (blauw).

aan ons voorbijgaat.

2.2.2

Returns

In figuur 2.8 zijn de genormaliseerde logaritmische returns weergegeven voor de Duitse DAX index (groen), Gaussisch verdeelde returns (rood) en returns verdeeld volgens een LAS distributie (blauw). We zien direct dat de returns uit een GBB simulatie (rood) zich braaf gedragen ten opzichte van de reële DAX index (groen). Bij de returns van de DAX index komen grote uitschieters voor die grote verliezen (of winsten) voorstellen. De returns van de GBB simulatie met LAS randomvariabele (blauw) vertonen wel hoge pieken maar, en dit is een ander belangrijk punt dat moeilijk te simuleren is, de pieken zijn niet geclusterd zoals bij de DAX index.


16

Figuur 2.8: Een vergelijking van de genormaliseerde returns van de DAX index en een GBB. De DAX index (groen) samen met een GBB simulatie met de random variabele getrokken uit een gaussische distributie (rood) en uit een Lévy alpha stable distributie (µ = 0, 5, a = 1, β = 0, m = 0) (blauw)

De reden dat we ook gebruik gemaakt hebben van de Lévy alpha-stable distributie is te zien in de figuren 2.9 en 2.10. De reden dat hier ook de Lévy alpha-stable distributie gebruikt is, is dat de PDF in zijn staarten een machtswet volgt. We zien dit fenomeen ook terug bij de PDF van de DAX returns. Dit is het gemakkelijkst voor te stellen door een log-log plot van de cumulatieve distributiefunctie van de returns aangezien deze meer data bezit in de staarten en dus minder afhangt van de bin-grootte [12]. In figuur 2.10 zijn de beide CDF’s gefit aan een machtswet. In [12] wordt ook aan getoond dat als de PDF een machtswet volgt met exponent α, de CDF ook een machtswet volgt met exponent (α − 1). De fit aan de positieve staart geeft een machtswet

met exponent −3, 73 en de fit aan de negatieve staart geeft een machtswet met exponent


17

Figuur 2.9: De probabiliteitsdistributiefuncties (PDF’s) voor de slotkoersen van de DAX index (groen) in de periode 1991-2011, een GBB met gaussische randomvariabele (rood) en een GBB met LAS random variabele (µ = 0, 5, a = 1, β = 0, m = 0) (blauw)

−5, 94. We kunnen dus besluiten dat de PDF van de returns van de DAX index aan het

positieve uiteinde een machtswet volgt met exponent −4, 73 en aan het negatieve uiteinde een machtswet met exponent −6, 94. Dit betekent dat extreme events niet zo onwaarschijnlijk zijn

als in het gewone GBB model wordt aangenomen waar de PDF van de returns exponentieel daalt daar deze een Gaussische distributie is. Aangezien deze ’extreme events’ een belangrijke impact hebben op het gedrag en de posities van de beleggers, is het duidelijk dat het GBB model slechts een serieuze benadering is van de dynamica die de financiële markten drijft.

2.2.3

Volatiliteit

Vervolgens bekijken we de volatiliteit. In figuur 2.11 wordt de volatiliteit getoond van de GBB simulatie (rood) en de DAX index (groen). Om goed te kunnen vergelijken zijn beide genormaliseerd.


18

Figuur 2.10: De cumulatieve distributiefuncties (CDF) voor de beide helften ([−10, 0[ en [0, 10]) van de PDF (figuur 2.9) van de DAX index. De staarten van de CDF’s zijn gefit aan een machtswet y = Cxα . De exponent van de fit aan de positieve returns is −3, 73, die van de fit aan de negatieve −5, 94.

Weer is er een duidelijk verschil. De volatiliteit van de GBB simulatie schommelt weinig rond een bepaalde waarde wat te verwachten is aangezien de σ uit de SDV van het Îto-proces constant is. Bij de volatiliteit van de DAX index liggen de zaken anders. Tijdens sommige beursperiodes zijn er pieken in de volatiliteit te zien die wijzen op het feit dat de markt zeer onrustig is. Het moment waar de pieken zich voordoen doen er zich grote prijsstijgingen en -dalingen voor, te zien op figuur 2.8. Het chaotisch karakter van de volatiliteit van de DAX index wijst erop dat de markt grotendeels niet in evenwicht is, in tegenstelling wat er in de EMH beweerd wordt. Het eenvoudige GBB model kan de complexe dynamica van de reële markten niet goed nabootsen.


19

Figuur 2.11: Een vergelijking van de genormaliseerde volatiliteit van de GGB (rood) en van de DAX (groen) index.

2.2.4

Informatie¨ entropie

Entropie zal later nog een belangrijk begrip worden dus laten we ook eens kijken naar de evolutie van de informatieëntropie. Deze is bepaald door H(τ ) =

X b

−pb [(τ − t) → τ ] ln pb [(τ − t) → τ ] ,

(2.24)

De som gaat over alle bins. pb is kans dat de return zich in de b’de bin bevindt. Met pb [(τ − t) → τ ] wordt bedoeld dat de discrete PDF berekend is uit het histogram van de

returns S(τ − t) tot S(τ ). Om een relevante waarde uit te komen voor de entropie moeten

er genoeg returns worden opgenomen om een goede PDF te vormen. Om de entropie op een

bepaald tijdstip te berekenen zullen we telkens de PDF van 200 returns gebruiken. Het is dus de informatieëntropie van de slotkoersen van 200 beursdagen.


20

Zoals eerder gezegd, is het belangrijk een goede bingrootte te kiezen. Om de entropieën van de GBB en de DAX index te vergelijken is er gekozen om bij het maken van het histogram het interval [−10, 10] te verdelen in 1000 bins met elk een grootte van 0, 02 In figuur 2.12 is het verloop van de informatieëntropie van de GBB simulatie en de DAX index gegeven. Wat opvalt is dat de dynamica van de entropie zeer gelijkaardig is aan die van de volatiliteit. In figuur 2.13 worden beiden boven elkaar geplot en is de overeenkomst duidelijk te zien. Deze gelijkenis geldt echter niet voor elke keuze van T . Voor T = 30 bv. vinden we geen goede overeenkomst tussen de entropie en de volatiliteit. Figuur 2.12: Vergelijking genormaliseerde Informatieëntropie GBB (rood) en DAX index (groen). De PDF’s van de genormaliseerde returns werden berekend in het interval [−10, 10] met bingrootte 0, 02.


21

Figuur 2.13: Entropie (groen)en volatiliteit (rood) van de DAX index in de periode 1991-2011. De dynamica van de entropie is zeer gelijkaardig aan die van de volatiliteit.

De entropie kan nog andere informatie geven. Als we entropie berekenen met slechts 2 bins, één van [−10, 0[ en één van [0, 10] dan krijgen we figuur 2.14. In sectie 2.1.9 hebben we gezien dat de entropie van een Gaussische distributie maximaal is. De afwijking in de entropie van de DAX index is er ten gevolge van een afwijking in de symmetrie van de PDF van de returns van de DAX index. Deze daling van entropie geeft dus weer dat de marktindexen aan het stijgen of dalen zijn met andere woorden dat de markt niet in evenwicht is. Dit komt overeen met de interpretatie van de entropie in de thermodynamica: een systeem dat niet in evenwicht is, heeft een lagere entropie dan een systeem in evenwicht, waarvoor de entropie steeds maximaal is.


22

Figuur 2.14: De genormaliseerde informatieëntropie van de DAX index (groen) berekend met 2 bins vergeleken met de entropie van een Gaussische distributie. Zoals voorspeld is de entropie van de Gaussische distributie steeds maximaal.

2.2.5

Autocorrelatie

De autocorrelatie van de returns van GBB simulatie en de DAX index zijn beiden gelijkaardig. De returns van de sluitingsprijzen voor de beursdagen zijn zo goed als niet gecorreleerd. Een hoge return op het ene moment verraadt niet wat er het volgende moment gaat gebeuren. Dit zit impliciet in het GBB model daar het Wiener proces een Markov proces is. De correlatie van de volatiliteit is in beide gevallen wel verschillend. Uit de DAX index volgt dat de volatiliteit op een bepaald tijdstip nog een lange tijd invloed heeft op de toekomst. De correlatie van de volatiliteit van de GBB simulatie daalt veel sneller naar 0 en onderschat dus de tijd waarin de volatiliteit gecorreleerd is. We kunnen hieruit besluiten dat de markten verschillende tijdschalen hebben. Een zeer kleine voor de returns en een grote tijdschaal voor de volatiliteit.


23

Figuur 2.15: Autocorrelatiefunctie van de DAX index in de periode 1991-2011. De returns (blauw) zijn niet gecorreleerd aangezien de correlatiefunctie direct naar nul valt. De volatiliteit (rood) daarintegen vertoont wel een lange tijds correlatie die maar zeer traag daalt.

2.2.6

Conclusie

De vergelijking toont aan dat de dynamica van het financiële systeem veel ingewikkelder in elkaar zit dan het eenvoudige GBB model. Er kan geopperd worden dat de uitschieters van de returns het gevolg zijn van gebeurtenissen buiten de financiële wereld (natuurrampen, terreuraanslagen, ...) maar dit is zeker niet altijd het geval (bv. [11]). Ook de clustering van extreme returns en hoge volatiliteit kan niet bekomen worden uit een GBB simulatie. Deze clustering wijst op een complexere dynamica dan dat er in de vergelijking van het Itô-proces gevat zit. De machtswetstaarten van de probabiliteitsdistributiefunctie kunnen wijzen op het bestaan van kritische fenomenen in de financiële markt (self-organized criticality? een soort fase transitie?). In de fysica zijn deze fenomenen met een machtswetstaarten distributie al uitvoerig onderzocht (zie [12], [13]). Er kan dus besloten worden dat het voorkomen van


24

Figuur 2.16: Autocorrelatiefunctie van de GBB simulatie. De returns (blauw) zijn niet gecorreleerd. Deze eigenschap zit impliciet in het stochastisch proces. Het Wiener proces is immers een Markov proces. De volatiliteit vertoont een grotere correlatie maar valt redelijk snel lineair naar nul.

crashes in de dynamica van het systeem zit en het is een uitdaging om een model te vinden dat dit zo goed mogelijk beschrijft. Een andere interessante conclusie is dat de volatiliteit en de entropie zich gelijkaardig gedragen en we dus in plaats van de volatiliteit ook de entropie kunnen bekijken om de onrust op de markten te onderzoeken.

Hoofdstuk 3

Moleculaire Dynamica Om de dynamica van de financiële markten te simuleren, gaan we nu gebruik maken van een Moleculaire Dynamica (MD) simulatie. De dynamica van de fysische grootheden uit een fysisch systeem zoals temperatuur, energie en entropie, kunnen bepaald worden door de microscopische bewegingsvergelijkingen numeriek te integreren. F~i (r1 , r2 , · · · , rN ) d2 r~i (t) = dt2 m

(3.1)

Eens de posities en de snelheden van elk deeltje in de simulatie gekend is, kunnen de macroscopische grootheden berekend worden. Merk op dat dit een klassieke benadering is. De bewegingsvergelijkingen zijn die van Newton en de invloed van de deeltjes op elkaar wordt voorgesteld door een empirische potentiaal. Om deze klassieke benadering te laten gelden, worden als deeltjes vaak argon atomen genomen. Kwantumeffecten kunnen verwaarloosd worden wanneer [14] λdB << l , 1

V 3 ) de gemiddelde afstand is tussen de deeltjes en λ = waarbij l = ( N

(3.2) √ h 3mkT

de de Broglie

golflengte. Door de hoge massa van het argon atoom en bij een temperatuur van 300K is de de Broglie golflengte voor een argon atoom van de orde 10−7 ˚ A. We zullen bij elke simulatie het volume V van het systeem zo aanpassen dat de gemiddelde afstand tussen de atomen ˚ is. Hierdoor is er altijd ruim aan de voorwaarde voldaan. Ook moet er geen van de orde 1A rekening gehouden worden met excitaties van het atoom. De excitatieënergie van een argon atoom is van de orde 10eV , terwijl de gemiddelde kinetische energie bij kamertemperatuur van de orde 0.1eV is. Een klassieke benadering waarbij het argon atoom wordt voorgesteld als puntdeeltje en een bepaalde potentiaal U bezit, is dus gerechtvaardigd. Wegens beperkte computercapaciteit kunnen we enkel systemen met een laag aantal deeltjes simuleren. Dit aantal 25

26

Hoofdstuk 3. Moleculaire Dynamica

deeltjes is 20 grootordes kleiner dan een systeem zoals ze normaal in de natuur voorkomen (orde 1023 deeltjes). Om randeffecten te vermijden zullen we gebruik maken van periodieke randvoorwaarden waardoor het bestudeerde systeem wordt omringd door een groot aantal kopies van zichzelf.

3.1

Oplossingsmethode

In een MD simulatie worden op elk tijdstip de bewegingsvergelijking ge¨ıntegreerd en zo de posities en de snelheden van elk deeltje op elk tijdstip bepaald. Hiervoor bestaan verschillende algoritmes, maar in dit werk zullen we gebruik maken van het Velocity Verlet algoritme. Verder hebben we ook nog een uitdrukking nodig voor de interactiepotentiaal U die de krachten tussen de deeltjes vastlegt.

3.1.1

Verlet algoritme

De meest gebruikte methode die gebruikt wordt om de Newtoniaanse bewegingsvergelijkingen op te lossen, is het Velocity Verlet algoritme. Met dit algoritme kunnen de posities en de snelheden van de deeltjes op elke tijdstap berekend worden uit hun vorige waarden. 1 r~i (t + ∆t) = r~i (t) + v~i (t)∆t + a~i ∆t2 , 2 a~i (t) + a~i (t + ∆t) v~i (t + ∆t) = v~i (t) + ∆t . 2

(3.3) (3.4) (3.5)

waarbij de versnellingen a~i (t) kunnen bepaald worden uit F~i (t) . (3.6) mi Een groot voordeel van dit algoritme is dat het zeer lang stabiel is (energie wordt behouden) a~i (t) =

waardoor er lang gesimuleerd kan worden.

3.1.2

Interactiepotentiaal

Het gedrag van het systeem in afhankelijk van de gekozen interactiepotentiaal. Hieruit worden de krachten berekend die de deeltje op elkaar uit oefenen. Een veel gebruikte empirische potentiaal waaruit de krachten kunnen afgeleid worden is de Lennard-Jones (LJ) Potentiaal. ULJ (r) = 4

σ 12 r

−

σ 6 r

.

(3.7)

27


met de diepte van de potentiaalput en σ de afstand waarop de potentiaal 0 wordt. Beide zijn afhankelijk van het soort deeltjes in het systeem. De r−12 -term stelt de elektrostatische repulsie voor en de r−6 -term is het gevolg van de aantrekkend Van Der Waals interacties. De Lennard-Jones potentiaal is een zeer goede benadering van de experimenteel waargenomen potentiaal tussen twee atomen of moleculen. Deze potentiaal is zeer goed geschikt voor MD van systemen in evenwicht. Als we, zoals hier in deze thesis, het systeem uit evenwicht brengen kan dat resulteren in deeltjes die een zeer grote afstand ∆x afleggen. Hierdoor kunnen deeltjes elkaar zeer dicht naderen en kan zo de potentiële energie zeer groot worden wat tot een onstabiele simulatie kan leiden [15]. Een goed alternatief voor de Lennard-Jones potentiaal die deze problemen niet heeft, is de Softcore (SC) potentiaal [16]. USC (r) =

(r − RA )2 H . − UA exp − 2 1 + exp ∆(r − Rr ) 2δA

(3.8)

De SC potentiaal behoudt de lange dracht interacties van de LJ potentiaal, maar zorgt voor een afvlakking van het centrale gedeelte, de ’harde kern’. Vandaar de naam ’Softcore’. De parameters ∆, δA , RA en Rr worden zo gekozen dat het lange dracht gedeelte overeenkomt met de LJ voor argon atomen. H wordt de hardheidsparameter genoemd en bepaalt de hoogte van de potentiaalbarrière. De waarden van de parameters worden weergegeven in tabel 3.1 en in figuur 3.1 worden de twee verschillende potentialen vergeleken.

3.1.3

Systeemheden

Om de berekeningen te vereenvoudigen werkt men met gereduceerde eenheden. De eenheden worden gereduceerd door ze te delen door een karakteristieke lengte- en energieschaal. Deze schalen halen we uit de LJ potentiaal. σ 12 σ 6 − ULJ (r) = 4 r r

(3.9)

Hier is de karakteristieke energieschaal en σ de karakteristieke lengteschaal. Voor Argon atomen worden de parameters ∆

2 2δA

RA

Rr

H

UA

39,4

0,062

1

0,79

20

1

Tabel 3.1: Parameters (in systeemeenheden) van gefitte SC potentiaal [15]


28

Figuur 3.1: Een vergelijking van de Lennard-Jones potentiaal en de Softcore potentiaal met parameters uit tabel 3.1. De afstand r staat uitgedrukt in systeemeenheden.

mAr = 6, 633517 · 10−26 kg Ar = 119, 8K kB σA = 3, 405 · 10−10 m

De gereduceerde eenheden worden dan

29


E Ar r ∗ r = σAr r Ar ∗ t = 2 t mAr σAr kB T∗ = T Ar

E∗ =

Alle figuren en grafieken over de gegevens van de MD simulatie zullen uitgedrukt zijn in deze gereduceerde eenheden.

3.2 3.2.1

Programma Initialisatie

Eerst en vooral worden het aantal deeltjes N, de dichtheid ρ, de gewenste temperatuur T en het aantal tijdstappen van de simulatie ingegeven. De startposities van de Argon atomen liggen op een FCC rooster daar deze configuratie het meest stabiel is. Vervolgens wordt aan elk deeltje een snelheid gegeven zodat de snelheden verdeeld zijn volgens een Maxwell-Boltzmann distributie. We doen dit door de snelheden te trekken uit een Gaussische distributie met q kB T gemiddelde 0 en standaardafwijking m . Deze snelheden worden nog gecorrigeerd door de totale impuls van het systeem te berekenen en van elke snelheid

ptot ~ Nm

af te trekken zodat de

totale impuls van het systeem nul is.

3.2.2

Thermalisatie

Na de initialisatie worden de bewegingsvergelijkingen ge¨ıntegreerd met behulp van het Verlet algoritme. De temperatuur van het systeem zal evenwel nog afwijken van de gewenste evenwichtstemperatuur Tev . Om dit recht te zetten wordt na een bepaald aantal tijdstappen (in ons geval 20) de snelheden herschaald met een factor λ.

λ= Volgens het equipartitietheorema

1

s

3kB (N − 1)Tev P . 2 i mvi

(3.10)

zal als λ naar 1 convergeert, T naar Tev convergeren.

Als deze waarde voor λ met voldoende nauwkeurigheid bereikt is en de fluctuaties op de 1

Het equipartitietheorema zegt:

P

i

mvi2 = 3kB (N − 1)T

30


opeenvolgende λ’s verwaarloosbaar worden dan is het evenwicht op de gewenste temperatuur bereikt en kan de simulatie van start gaan.

3.2.3

Productiefase

Tijdens de simulatie worden de belangrijke variabelen weggeschreven. Wat volgt is een collectie figuren met de resultaten van een MD simulatie met de SC potentiaal op een vaste temperatuur T ∗ = 1, 1. 3.2.3.1

Temperatuur

De temperatuur kan op elk moment berekend worden uit de snelheden van de deeltjes. Met behulp van het equipartitietheorema vindt men T (t) =

P

mvi2 (t) 3kB (N − 1)

(3.11)

31


Figuur 3.2: Het verloop van de temperatuur tijdens een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. Zoals verwacht blijft de temperatuur ongeveer dezelfde in dit gesloten systeem.

3.2.3.2

Energie

De kinetische en potentiële energie worden op elk tijdstip bepaald als

N

1X mvi2 (t) , Ek (t) = 2 Ep (t) =

i=1 N X

i<j=1

U (|~ ri − r~j |) .

(3.12)

(3.13) (3.14)

Als controle voor de stabiliteit van de MD simulatie bereken we ook de totale energie Etot (t) = Ek (t) + Ep (t). Deze moet constant zijn doorheen de simulatie.

32


Figuur 3.3: De energieën van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5.

3.2.3.3

∆x

∆ is de gemiddelde verandering van de positie van de deeltjes op opeenvolgende tijdstippen. Deze wordt berekend door voor elk deeltje de verandering ten opzichte van de vorige positie te bepalen en daarvan het gemiddelde te bepalen. N 1 X ∆x = ∆xi . N

(3.15)

i

Dit is de variabele die we gaan relateren aan de returns uit de marktdata. De returns uit de markt geven een prijsverandering weer, net zoals ∆x de positieverandering weergeeft. Om een relevante vergelijking te maken, zorgen we ervoor dat beiden genormeerd zijn. S(t) − hS(t)i ∆x(t) − h∆x(t)i ↔ σS σ∆x

(3.16)


33

Figuur 3.4: Het verloop van de gemiddelde verandering van de deeltjes van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5.

Uit deze variabele wordt, net zoals in Hoofdstuk 2, de entropie en de volatiliteit bepaald.

34


Figuur 3.5: Probabiliteitsdistributiefunctie ∆x van de deeltjes van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. Zoals verwacht is dit in evenwicht een gaussische distributie.

3.2.3.4

Entropie

Met entropie bedoelen we hier enkel de translatieëntropie (zie hoofdstuk Entropie). Dit doen we omdat de entropie dan op een gelijkaardige manier berekend wordt dan de informatieëntropie uit de financiële data. Aangezien de massa gewoon een evenredigheidsconstante is (die ook gewoon 1 kan verondersteld worden), kan de impuls gerelateerd worden aan ∆x (p ∝ ∆x). De translatieëntropie wordt in het discrete geval gegeven door Str = −kB

X

pi ln pi .

(3.17)

i

Hier is pi de kans dat een bepaalde ∆x in bin i is gelegen. In de praktijk zullen we weer het verloop van de entropie in de tijd berekenen.


Str (τ ) =

X i

−pi [(τ − T ) → τ ] ln pi [(τ − T ) → τ ] .

35

(3.18)

De waarde van T is niet zo gemakkelijk te kiezen aangezien dat 1 dag markttijd niet zomaar gelijk kan worden verondersteld als 1 simulatiestap. Het mappen van simulatietijd op markttijd zal later nog besproken worden. Voor een simulatie in evenwicht krijgen we de volgende figuur (hier T = 1.000). Figuur 3.6: Entropie in functie van de tijd van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. Zoals verwacht uit de thermodynamica schommelt de entropie een weinig rond een constante waarde voor een systeem in evenwicht.

36


3.2.3.5

Volatiliteit

Naar analogie met de financiële markt definiëren we ook een volatiliteit over een bepaald tijdstip als σt (τ ) =

Pτ

i=τ −T (S(i)

− µ)2

. (3.19) T −1 Om dezelfde reden als bij de entropie is ook hier de keuze van T niet triviaal. In figuur 3.7 is de volatiliteit weergegeven van een MD simulatie in evenwicht. Figuur 3.7: Volatiliteit in functie van de tijd van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5.

3.2.3.6

Snelheidsautocorrelatiefunctie

Om enig inzicht te krijgen in het mappen van markttijd op simulatietijd kunnen we gebruik maken van de snelheidsautocorrelatiefunctie. In hoofdstuk 2 hebben we gezien dat de autocorrelatiefunctie van de returns (figuur 2.15) na 1 dag al onmiddellijk naar nul valt. Met


37

andere woorden de returns zijn ongecorreleerd. Met behulp van de snelheidsautocorrelatiefunctie kunnen we bepalen na hoeveel tijdstappen de snelheden (∆x’s) ongecorreleerd zijn. We kunnen dan stellen dat 1 dag markttijd overeenkomt met het aantal tijdstappen dat nodig is om de snelheden ongecorreleerd te maken. In figuur 3.8 is de snelheidscorrelatiefunctie gegeven van een MD simulatie in evenwicht. Figuur 3.8: Snelheidsautocorrelatiefunctie van een MD simulatie met 256 deeltjes en een initiëele temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5.

We zien dat de ∆x’s ongecorreleerd zijn na afgerond 300 tijdstappen. We kunnen dus ruwweg stellen dat 1 marktdag overeenkomt met 300 simulatiestappen.

Hoofdstuk 4

Entropie Entropie. Een zeer gekende grootheid die in veel ingenieurs- en wetenschapsdisciplines een betekenis heeft maar waarvoor geen eenduidige betekenis bestaat. Entropie werd voor het eerste fenomelogisch gedefinieerd in de vorm van de 2de hoofdwet van de thermodynamica die zegt dat de entropie van een ge¨ısoleerd systeem enkel kan stijgen of constant blijft. Het gaf een grootheid aan een feit dat door vele wetenschappers werd vastgesteld dat in ge¨ısoleerde systemen de spontane processen zich in 1 richting begeven, denk maar aan het Carnot-proces. Bij de analyse van de financiële markt hebben we gebruik gemaakt van de informatiëentropie gedefinieerd door Shannon zogeheten door de intu¨ıtieve link met het entropie begrip in de statistische mechanica. Zoals Shannon het ooit zei: ” My greatest concern was what to call it. I thought of calling it information, but the word was overly used, so I decided to call it uncertainty. When I discussed it with John von Neumann, he had a better idea. Von Neumann told me, You should call it entropy, for two reasons. In the first place your uncertainty function has been used in statistical mechanics under that name, so it already has a name. In the second place, and more important, nobody knows what entropy really is, so in a debate you will always have the advantage. ” Maar in hoeverre kunnen we de informatieëntropie van een markt vergelijken met de statistische entropie bepaald uit een fysische simulatie? Kan er hoegenaamd sprake zijn van een consistente vergelijking? Op deze vragen zullen we uiteindelijk een antwoord willen vinden.

38

39

Hoofdstuk 4. Entropie

4.1

Entropie in de Statistische Mechanica

De entropie als grootheid voor het statistisch beschrijven van systemen, werd ingevoerd door Ludwig Boltzmann in de bekende formule S = kb ln Ω .

(4.1)

Hier is Ω het aantal microtoestanden waarin een systeem zich kan bevinden. Deze formule is enkel geldig als we over systemen in evenwicht spreken. Een meer algemene formule werd opgesteld door J.W. Gibbs S = −kb

X

Pj ln Pj .

(4.2)

j

Pj is de probabiliteit dat het systeem zich in een bepaalde microtoestand j bevindt. Als alle microtoestanden even waarschijnlijk zijn (Pj =

1 Ω)

dan herleidt 4.2 zich tot 4.1.

Om Pj te bepalen beschouwen het canonisch ensemble 4.1. Het canonisch ensemble wordt gedefinieerd als een ensemble dat is samengesteld uit een aftelbaar aantal subsystemen die van elkaar gescheiden zijn door een vaste wand die ondoorlaatbaar is voor deeltjes maar wel warmte doorlaat. Elk subsysteem bevindt zich dus in een warmtebad gevormd door de andere subsystemen. Voor elk systeem zijn de temperatuur T, het volume V en het aantal deeltjes N constant. Tussen de verschillende subsystemen kan energie worden uitgewisseld. De energie van het gehele ge¨ısoleerde systeem is echter een constante Etot . We kijken nu naar een subsysteem j dat in contact staat met een reservoir (de rest van de subsystemen in het ensemble). We zoeken nu naar de probabiliteit Pj dat dit subsysteem een energie heeft van Ej . Aangezien het subsysteem een energie heeft van Ej kan het reservoir zich in elke microtoestand bevinden met een energie Er = Etot − Ej . Volgens het ergodetisch principe zijn elk van deze microtoestanden even waarschijnlijk. Hieruit volgt dat de probabi-

liteit dat het subsysteem j een energie Ej heeft evenredig is met het aantal microtoestanden Ωr (Etot − Ej ) die het reservoir kan aannemen. Met behulp van de Boltzmannentropie 4.1 kunnen we dit omvormen tot

Pj ∝ Ωr (Etot − Ej ) , Sr (Etot − Ej ) . ∝ exp kb

(4.3) (4.4)

met Sr (Etot − Ej ) de entropie van het reservoir. We veronderstellen dat de energie van het subsysteem j zeer klein is ten opzichte van de totale energie van het systeem wat een goede

40


Figuur 4.1: Het kanonisch ensemble (figuur uit [17])

benadering is zolang het systeem in evenwicht is en de energiefluctuaties niet te groot zijn. Met Ej <<< Etot kunnen we de Taylor ontwikkeling van Sr (Etot − Ej ) rond Etot afbreken na de 2de term.

Sr (Etot − Ej ) ≈ Sr (Etot ) − Ej Nu is

∂Sr (Etot ) ∂Etot

=

1 T

∂Sr (Etot ) . ∂Etot

(4.5)

met T de temperatuur van het reservoir. Sr (Etot ) is een constante zodat

we vinden (met β =

1 kb T )

Ej Pj ∝ exp Sr (Etot ) − kb T

,

(4.6)

∝ exp (−βEj ) . De proportionaliteitsconstante vinden we door op te merken dat

(4.7) P

j

Pj = 1. Hieruit volgt

41


Pj = C exp (−βEj ) ,

Z=

P

4.2

j

(4.8)

exp (−βEj ) =P . j exp (−βEj )

(4.9)

exp (−βEj ) wordt de Zustandssumme of de partitiefunctie genoemd.

Entropie van veeldeeltjessysteem

Met behulp van de uitdrukking voor Pj kunnen we nu de entropie van een systeem met een vast aantal deeltjes en een vaste temperatuur in contact met een warmtebad bepalen. We zullen bespreken hoe de entropie van een N deeltjessysteem wordt berekend in de statistische mechanica. De Hamiltoniaan van N identieke deeltjes is

H=

N X X p2i + U (i, j) . 2m

(4.10)

i<j

i

Eerst bepalen we de partitiefunctie ZN van de N deeltjes. ZN =

X

exp (−βH) .

(4.11)

De som loopt over alle mogelijke microtoestanden. Elke microtoestand wordt gedefinieerd door een bepaald punt in de faseruimte die zoals we weten uit de kwantummechanica, gra1 V 3 kunnen we de nulair is. We werken echter klassiek en als er geldt dat λdB <<< N

thermodynamische limiet toepassen en de sommatie vervangen door een integraal. ZN met K =

P

p2i i 2m

1 = 1 N !h3N

Z

d3N pd3N r exp (−β(K + U )) .

de kinetische energie term en U =

P

i<j

(4.12)

Ui,j de potentiële energie term.

De partitiefunctie valt nu uiteen in 2 delen, 1 met de kinetische energie term en een met de potentiële energie term. tr int ZN = ZN × ZN .

(4.13)

Beide termen zullen apart behandeld worden en elk hun bijdrage leveren aan de entropie.

4.2.1

Translatie¨ entropie

De translatieëntropie is de entropie van een systeem waar de deeltjes onderling niet interageren. Dit is de definitie van een ideaal gas. De energie van het systeem bestaat enkel uit

42


kinetische energie, enkel afhankelijk van de impuls van de deeltjes. Klassiek is de uitdrukking voor de entropie S

tr

= −kb

Z

P tr (~ p) ln P tr (~ p) .

P exp −β N i

P tr =

pi 2 2m

(4.14)

. (4.15) Z De deeltjes zijn niet onderscheidbaar waardoor de berekening van de partitiefunctie een pak eenvoudiger wordt 1 tr [Z (T, V )]N . (4.16) N! 1 De één-deeltjespartitie-functie van een niet-interagerend deeltje kan eenvoudig berekend wortr ZN (T, V, N ) =

den [18]

Z1tr (T, V

)=

Z

∞

0

=

V . λ3

p2 4πp2 dp exp −β , h3 2m

(4.17) (4.18)

met λ de thermische golflengte λ= √

h . 2πmkT

(4.19)

Hierdoor wordt Pjtr gegeven door N

Pjtr

X pi 2 N !λ3N = exp −β VN 2m i

!

,

(4.20)

zodat de translatieëntropie kan berekend worden

VN S tr (~ p) = −kb 3N h N!

Z

N

X pi 2 N !λ3N d3 p~1 d3 p~2 · · · d3 p~N exp −β VN 2m i

!

"

N

X pi 2 N !λ3N ln exp −β VN 2m i

!" # 3N Z N N X X λ pi 2 N !λ3N pi 2 3 3 3 = −kb d p~1 d p~2 · · · d p~N exp −β ln −β . h 2m VN 2m i

(4.21)

i

(4.22)

Dit zijn N onafhankelijke integralen en kunnen dus gemakkelijk1 berekend worden. We vinden 1

Met behulp van de bekende Gaussische integralen In (a) =

R∞ 0

dxxn exp(−ax2 ) met I2 (a) =

1 4a

π a

1 2

.

!#

,

43


S tr = −kb ln

N !λ3N VN

+ kb

3N . 2

(4.23)

Beschouwen we een systeem met een groot aantal deeltjes dan vereenvoudigt S tr zich met behulp van de stirling formule2 nog tot S

4.2.2

5 V + . = N kb ln N λ3 2

tr

(4.24)

Interactie¨ entropie

Als het systeem bestaat uit interagerende deeltjes waarbij de interactie tussen twee deeltjes beschreven wordt door de potentiaal U (i, j), krijgen we voor de probabiliteitsdistributie en de toestandssom de volgende uitdrukkingen

P int =

Z int = =

1 h3N N ! V

P exp −β i<j U (i, j) Z int

Z

N

h3N N !

ζ.



d3N r exp −β

X i<j

,

(4.25) 

U (i, j) ,

(4.26) (4.27)

met 1 ζ= N V De interactieëntropie gegeven door S

int

Z



d3N r exp −β

= −kb

Z

X i<j



U (i, j)

(4.28)

P int (~r) ln P int (~r) ,

(4.29)

wordt, na invulling van de PDF (4.25) en de partitiefunctie (4.27)

S int (~r) = −kb

Z

d3N r

1

h3N N !

h3N N !

VN

V Nζ = kb ln 3N − h N! 2

ln n! ≈ n ln n − n

Z

P exp −β i<j U (i, j) ζ



ln 

h3N N ! exp VN

P exp −β i<j U (i, j) X 1 d3N N [−β U (i, j)] . V ζ i<j

−β

 U (i, j) i<j , ζ

P

(4.30)

(4.31)

44


De tweede term in deze uitdrukking kunnen we nog herschrijven.

P

i<j

U (i, j) is niet anders

het gemiddelde van de potentiële energie U = hEp (t)i. De integraal is gelijk aan 1 ten gevolge van de normalisatieconditie:

Z

d3N r

1

P (r) , P Z exp −β U (i, j) i<j 1 = d3N N . V ζ

1=

h3N N !

(4.32) (4.33) (4.34)

Hierdoor vereenvoudigt de interactieëntropie zich tot V Nζ + kb βU . (4.35) h3N N ! U kan berekend worden uit de potentiële energie. Voor ζ is echter geen gesloten uitdrukking S int = kb ln

te vinden. Hierom zullen we ζ berekenen met behulp van een perturbatie-expansie. Hiervoor worden de hulpfuncties van Mayer [19] voorgesteld: λ(i, j) = exp (−βU (i, j)) − 1 .

(4.36)

Deze worden gebruikt om de volgende reeksontwikkeling op te stellen. Y

(1 + λ(i, j)) = 1 +

i<j

X

λ(k, l) +

k
X

k
λ(k, l)λ(n, m) + · · ·

(4.37)

Deze reeksontwikkeling is gemakkelijk te interpreteren. De eerste term geeft de ideale gas bijdrage, terwijl de tweede term de interactie tussen twee nabij gelegen deeltjes weergeeft. De derde term geeft de interactie weer tussen drie nabij gelegen deeltjes. Deze term zal maar belangrijk worden bij hogere dichtheden. Aangezien we de MD simulatie enkel zullen uitvoeren in de gas- en vloeistof fase, kunnen we de reeksontwikkeling na de tweede term afbreken. 1 ζ= N V

Z

d3N r

1 ≈1+ N V

Z

Y

exp (−βU (i, j)) ,

(4.38)

X

(4.39)

i<j

d3N r

λ(k, l) .

k
De berekening van de integraal is sterk vereenvoudigd daar we bij elke term van de som slechts rekening moeten houden met 2 plaatsvectoren. Met deze uitdrukking kunnen we dus de interactieëntropie van N interagerende deeltjes berekenen.

45


4.3

Entropie in de Simulatie

In hoofdstuk 2 is de informatieëntropie van marktdata berekend uit de PDF van de returns. Mappen we nu de gemiddelde afstand die de deeltjes tijdens de simulatie afleggen (∆x) op de returns, dan kunnen we de entropie in de simulatie op een vergelijkbare manier berekenen als bij de marktdata. De entropie van de simulatie wordt dan berekend uit de PDF van de stapgroottes ∆x zoals gedaan is in vergelijking 3.18. Is deze intu¨ıtieve manier om de entropie uit de MD simulatie te berekenen, gevest op enige fysische grond? Zoals we gezien hebben wordt S tr berekend met de kans Pj dat het systeem zich in een bepaalde energietoestand Ej bevindt. Met de gelijkheid E =

p2 2m

kunnen we de kans Pj relateren aan de PDF van de

stapgrootte daar m in de simulatie gelijk is aan 1. De kans Pj kunnen we schrijven als

Pjtr

X p2 i ∝ exp β 2m i 2 vtot . ∝ exp 2kb T

!

,

(4.40) (4.41) (4.42)

2 = Hierin is vtot

P

2 i vi .

Voor een systeem in evenwicht is de PDF van de stapgroottes een

gaussische distributie (zie figuur 3.5). Pj∆x

∆x2 ∝ exp − 2 2σ

(4.43)

De entropie die berekend is uit de PDF van de stapgroottes komt dus overeen met de translatieëntropie van het systeem. We zullen dus in deze thesis enkel gebruik maken van de translatiëntropie van het fysische systeem om een vergelijking te maken met de informatieëntropie van de financiële markten. Hieruit volgt ook dat de vorm van de potentiaal U geen rechtstreekse invloed heeft op de waarde van de entropie die vergeleken wordt met de informatieëntropie uit de marktdata aangezien we de aanwezigheid van de interactieëntropie niet in rekening brengen.

Hoofdstuk 5

Simuleren van Markten met NEMD In dit hoofdstuk zullen we proberen met behulp van een MD simulatie die uit evenwicht wordt gebracht de dynamica van de markt te simuleren. Het moleculair systeem zal uit evenwicht worden gebracht door met de temperatuur te spelen. Eerst zullen we onderzoeken hoe het systeem reageert op plotse en geleidelijke temperatuursveranderingen. Met deze inzichten kunnen we dan proberen een simulatie uit te voeren die dezelfde dynamica bevat als een reële markt.

5.1

Temperatuurschok

Eerst gaan we eens kijken wat er gebeurt als we het systeem plots op een andere temperatuur brengen. Dit wordt in de praktijk gedaan door op een bepaalde tijdstap de snelheden van alle deeltjes te vermenigvuldigen met λ.

λ=

s

3kB (N − 1)Tf P , 2 i mvi

(5.1)

met Tf gelijk aan de nieuwe temperatuur die we het systeem willen opleggen. De stijging in de volatiliteit en de entropie wordt veroorzaakt door de overgang van de ene evenwichtstoestand naar de andere. Tijdens de temperatuurschok veranderen de snelheden naar een andere waarde, waardoor het histogram van de plaatsverandering ∆x(t) geen gewone gaussische distributie is. In figuur 5.1 wordt het temperatuursverloop van een MD simulatie weergegeven waarvan op tijdstappen 20.000, 40.000, 60.000 en 80.000 de temperatuur plots sterk wordt verhoogd of verlaagd en dan op een nieuwe temperatuur wordt gezet. De invloeden van de temperatuursveranderingen worden weergegeven het histogram van ∆x(t) in figuur 5.2. De twee grote (gaussische) pieken komen van het systeem in evenwicht op de twee verschillende temperaturen. De grote temperatuurschokken geven aanleiding tot de staarten 46

Hoofdstuk 5. Simuleren van Markten met NEMD

47

Figuur 5.1: Het temperatuursverloop van een NEMD simulatie waarbij op tijdstappen 20.000, 40.000, 60.000 en 80.000 de temperatuur van het systeem plots wordt verhoogd of verlaagd.

van de probabiliteits distributie. Als ervoor kunnen zorgen dat deze staarten afvallen als een machtswet, dan kunnen we deze hoge en lage ∆x misschien linken aan de extreme events in de financiële markten.


48

Figuur 5.2: Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.1. De twee pieken stellen de ∆x’s voor in de twee evenwichtssituaties op temperatuur T ∗ = 0, 7 en T ∗ = 0, 95. Ten gevolge van de temperatuurschokken nemen we ook extremere waarden waar voor bepaalde ∆x.

Met deze PDF kunnen we nu ook de volatiliteit en de entropie berekenen. Zoals uitgelegd in het vorige hoofdstuk veronderstellen we dat 1 marktdag overeenkomt met 300 simulatiestappen. Als tijdsvenster voor het berekenen van de volatiliteit en entropie zullen we hier 10.000 tijdstappen gebruiken. In figuur 5.3 worden deze weergegeven. De stijgingen in de entropie en volatiliteit worden veroorzaakt door de plotselinge temperatuursverandering. Het is zo dat hoe drastischer de temperatuurswijziging, hoe sterker de verhoging in de volatiliteit en entropie is, wat intu¨ıtief duidelijk is. Een mogelijk manier om een reeks marktdata te simuleren is dus het systeem op verschillende tijdstappen een temperatuurschok te geven van verschillende groottes. Merk op dat de entropie en volatiliteit hetzelfde zijn als het systeem in evenwicht is, zelfs als de temperatuur anders is.


49

Figuur 5.3: De entropie en de volatiliteit in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.1. De stijging en daling in entropie en volatiliteit zijn het gevolg van de plotse temperatuursveranderingen.

5.2

Temperatuursgradi¨ ent

Nu gaan we het systeem een temperatuursgradiënt opleggen. Vanaf een bepaalde tijdstap wordt de temperatuur geleidelijk aan verhoogd. Aangezien we met een discrete simulatie werken (per definitie) kunnen we geen continue gradiënt opleggen. Ook is het niet mogelijk om de temperatuur op elk tijdstip te veranderen omdat we dan op elke tijdstap de temperatuur van het systeem vastleggen en dus de snelheden van de deeltjes. We kunnen dan niet meer zeggen dat dit een simulatie is van een systeem van deeltjes die met elkaar interageren, maar we zouden de deeltjes in een bepaalde toestand forceren gedurende de hele simulatie. Hierom verhogen we de temperatuur van het systeem pas om de 200 tijdstappen zodat het systeem zich toch even kan aanpassen aan de nieuwe temperatuur. De naam temperatuurgradiënt lijkt nu een beetje verkeerd gekozen. We onderwerpen het systeem eigenlijk aan een hele


50

reeks kleine temperatuurschokken met dezelfde groottte totdat het systeem zich op de nieuwe gewenste temperatuur bevindt. Figuur 5.4: Het temperatuursverloop van een NEMD simulatie waarbij tussen tijdstappen 10.000 en 20.000 de temperatuur van het systeem geleidelijk aan wordt verdubbeld.

In figuur 5.4 wordt het verloop van de temperatuur weergegeven. Doordat de temperatuurschokken tussen tijdstappen 10.000 en 20.000 klein genoeg zijn en toch voldoende dichtbijeen liggen zijn de temperatuurschokken niet zichtbaar op de grafiek, maar gaan verloren in de gaussische random ’noise’ van de MD simulatie. Het histogram van ∆x(t) in figuur 5.5 vertonen twee pieken afkomstig van de evenwichtssituaties voor en na de temperatuursverhoging. De invloed van de temperatuursverhoging zelf is te zien tussen de twee pieken. In plaats van de som van de twee nabij gelegen Gaussische pieken valt het histogram hier hoger uit door de temperatuursgradiënt.


51

Figuur 5.5: Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.4. De twee pieken stellen de ∆x’s voor in de twee evenwichtssituaties op temperatuur T ∗ = 0, 7 en T ∗ = 1, 3. Ten gevolge van de temperatuursgradiënt worden er waarden van ∆x’s waargenomen die zich tussen de twee evenwichtspieken bevinden.

De invloed op de entropie en de volatiliteit is duidelijk te zien in figuur 5.6. De lichte stijging is het gevolg van de temperatuurschokken waarvan de invloed zich 10.000 tijdstappen laat voelen.


52

Figuur 5.6: De entropie en de volatiliteit in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.4. De kleine stijging en daling in entropie en volatiliteit zijn het gevolg van de geleidelijke temperatuursveranderingen.

5.3

Marktsimulatie

Na deze twee experimenten om te zien hoe het systeem reageert op temperatuursveranderingen, zullen we proberen om een reële markt na te bootsen. Dit zullen we proberen door op bepaalde tijdstippen het systeem aan een temperatuurschok te onderwerpen. De grootte van de temperatuurschokken en de tijdstippen waarop deze plaatsvinden zijn nog onbekende parameters die we nog kunnen ’tunen’.

5.3.1

’Gaussische’ schokken

Eerst gaan we het volgende proberen. We veranderen de temperatuur van het systeem na elke 2.000 tijdstappen. 200 tijdstappen na het veranderen van de temperatuur zetten we de


53

temperatuur terug op zijn evenwichtswaarde. De grootte van de temperatuursverandering wordt bepaald door een gaussische randomvariabele. Concreet wordt de temperatuur na elke 2000 tijdstappen vermenigvuldigd met of gedeeld door (1 + |rnd|) (afhankelijk van het teken van de gaussische randomvariabele), met rnd een gaussische randomvariabele met gemiddelde

0 en standaardafwijking 1. Krijgen we als randomgetal 2, dan wordt de temperatuur met 3 vermenigvuldigd. Krijgen we daarentegen −1 als randomgetal dan wordt de temperatuur door 2 gedeeld. Er is gekozen om de temperatuur te vermenigvuldigen en te delen omdat

we dan nooit een negatieve temperatuur zouden kunnen uitkomen, wat wel het geval is als we zouden optellen en aftrekken. Het resultaat van wat deze veranderingen teweeg brengen worden getoond in figuur 5.7. Figuur 5.7: Het temperatuursverloop van een NEMD simulatie waarbij na elke 2.000 tijdstappen de temperatuur gedurende 200 tijdstappen wordt verhoogd of verlaagd met een grootte bepaald door een gaussische randomvariabele.

Het histogram van de snelheden in deze simulatie is weergegeven in figuur 5.8. Wat eerst


54

al opvalt is dat het histogram niet symmetrisch is. Dit komt door de implementatie van de temperatuursverandering. Een vermenigvuldiging met of deling door een factor 2 is beiden even waarschijnlijk. De absolute temperatuursverandering is echter niet gelijk. Beginnen we bijvoorbeeld met een temperatuur T ∗ = 1 dan krijgen we resp. 2 en 0, 5 als nieuwe temperaturen die een temperatuursverandering geven van in het ene geval 2 − 1 = 1 en in het

andere geval 0, 5 − 1 = −0, 5. Het is dus logisch dat het histogram niet symmetrisch is. Om de staarten te analyseren, bereken we de cumulatieve distributie functie (CDF). De positieve

staart wordt getoond in figuur 5.9, de negatieve in 5.10. Op beide figuren worden de CDF’s vergeleken met die van een gaussische distributie en met de DAX index (uit figuur 2.10). Figuur 5.8: Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.7. De asymmetrie in het histogram is een gevolg van de implementatie van de temperatuurschok.

De negatieve staart van de PDF van de ∆x’s volgt de staart van een Gaussische distributie. Dit is niet wat we willen vinden. De positieve staart volgt echter wel een interessant gedrag.


55

Figuur 5.9: De cumulatieve distributiefunctie van de positieve staart van het histogram in figuur 5.8 vergeleken met de CDF van een GBB simulatie en de CDF van de DAX index. Alhoewel grotere ∆x meer waarschijnlijk zijn, vertoond de positieve staart geen verval volgens een machtswet, maar vervalt eerder exponentieel.

Figuur 5.10: De cumulatieve distributiefunctie van de negatieve staart van het histogram in figuur 5.8 vergeleken met de CDF van een GBB simulatie en de CDF van de DAX index. De negatieve staart vertoont geen van de karakteristieken die we zoeken.

We zien duidelijk dat er zoals bij de DAX index zich extremen voordoen. Het is echter wel zo dat de staart afvalt zoals een exponentiële functie en niet zoals een machtswet. We zullen


56

dus iets anders moeten proberen. Voor de volledigheid wordt in figuur 5.11 de entropie en volatiliteit gegeven. Merk op dat deze veel minder chaotisch zijn als de volatiliteit en entropie uit de markten (figuren 2.11 en 2.12). De reden hiervoor zal later nog toegelicht worden. Figuur 5.11: De entropie en de volatiliteit in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.1. De stijging en daling in entropie en volatiliteit zijn het gevolg van de plotse temperatuursveranderingen.


5.3.2

57

’Levy’ schokken

We proberen nu hetzelfde trucje al in het vorige puntje. We gaan nu echt de randomvariabele die de grootte van de temperatuursverandering bepaald nu trekken uit een lévy alpha-stable (LAS) distributie met als parameters a = 0, 01, β = 0, µ = 0, 5 en m = 0. Hierdoor zullen we grotere temperatuursveranderingen induceren en misschien ook een vette staart in het histogram van de ∆x’s. Het temperatuursverloop wordt gegeven in figuur 5.13. Figuur 5.12: Het temperatuursverloop van een NEMD simulatie waarbij na elke 2.000 tijdstappen de temperatuur gedurende 200 tijdstappen wordt verhoogd of verlaagd met een grootte bepaald door een randomvariabele getrokken uit een LAS distributie met parameters a = 0, 01, β = 0, µ = 0, 5 en m = 0.

Om dezelfde reden als in de vorige paragraaf is het histogram van de snelheden niet symmetrisch. Het is al duidelijk uit de figuur 5.13 dat de positieve en negatieve staart sterk verschillen. Dit komt nog beter tot uiting in de figuren 5.14 en 5.15 van de CDF’s. De positieve staart (figuur 5.14) geeft een serieuze overschatting van de grote ∆x’s terwijl de kleine


58

∆x onderschat worden. Ook de negatieve staart (figuur 5.15)geeft geen goede overeenkomst met de CDF van de returns van de DAX index. Figuur 5.13: Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.7. De asymmetrie in het histogram is een gevolg van de implementatie van de temperatuurschok.

Het is duidelijk dat als we het systeem op deze manier uit evenwicht brengen we niet de gewilde effecten introduceren die de NEMD simulatie eenzelfde dynamica geeft als de markten. Wel is het zo dat bij deze simulatie de staarten van het histogram van de ∆x’s afloopt zoals een machtswet. We zoeken verder ...


59

Figuur 5.14: De cumulatieve distributiefunctie van de positieve staart van het histogram in figuur 5.8 vergeleken met de CDF van een GBB simulatie en de CDF van de DAX index. Alhoewel er een deel van de CDF vervalt als een machtswetswet, is dit niet het resultaat dat we zoeken.

Figuur 5.15: De cumulatieve distributiefunctie van de negatieve staart van het histogram in figuur 5.8 vergeleken met de CDF van een GBB simulatie en de CDF van de DAX index. Alhoewel de CDF vervalt als een machtswetswet, zijn de grote ∆x’s veel minder waarschijnlijk dan bij de CDF van de returns van de DAX index.


5.4

60

’Levy-Gauss’ schokken

In figuur 5.9 zien we dat we door het systeem te verstoren met temperatuurschokken bepaald door een gaussische randomvariabele een verloop van de CDF krijgen die zeer gelijkaardig is aan die van het verloop van de CDF van de DAX index. Er is echter 1 probleem: de staart valt niet af als een machtswet. Een mogelijk oplossing is niet alleen ’Gaussische’ schokken te gebruiken maar ook ’Levy’ schokken. In de volgende simulatie genereren we om de 1.000 tijdstappen een temperatuurschok die 70% kans heeft om een ’gaussische’ schok te zijn en 30% kans om een ’levy’ schok te zijn. Door het invoeren van een aantal ’levy’ schokken proberen we een machtswetstaart te genereren. Het temperatuursverloop is weergegeven in figuur 5.16. Figuur 5.16: Het temperatuursverloop van een NEMD simulatie waarbij na elke 1.000 tijdstappen de temperatuur gedurende 200 tijdstappen wordt verhoogd of verlaagd met een grootte bepaald door een randomvariabele voor 30% getrokken uit een LAS distributie met parameters a = 0, 01, β = 0, µ = 0, 5 en m = 0 en voor 70% getrokken uit een gaussische distributie met µ = 0 en σ = 1.


61

In het histogram in figuur 5.17 lijkt het erop dat de positieve staart een vette staart is. Dit is duidelijker te zien in figuur 5.18 waar de CDF van de positieve staart uit het histogram is weergegeven. We nemen een duidelijke machtswetstaart vast terwijl ook voor kleinere ∆x de CDF van de simulatie dicht aanleunt bij de CDF van de returns van de DAX index. Dit is wat we willen vinden en toont aan dat we met behulp van de juiste manipulaties de dynamica van een markt kunnen nabootsen in een NEMD simulatie. Figuur 5.17: Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.16. De asymmetrie in het histogram is een gevolg van de implementatie van de temperatuurschok.


62

Figuur 5.18: De cumulatieve distributiefunctie van de positieve staart van het histogram in figuur 5.8 vergeleken met de CDF van de DAX index. De CDF vertoont grote gelijkenissen met die van de DAX index.

5.5

Opmerkingen

Bij deze simulaties is het nodig nog een paar belangrijke opmerkingen te geven. Het enigszins mooie resultaat in figuur 5.18 is het beste resultaat uit een heleboel simulaties. Door het random karakter van de opgelegde temperatuurschokken en het feit dat de simulatie ’slechts’ 100.000 tijdstappen lang is, is het niet verwonderlijk dat er grote verschillen zijn tussen de verschillende simulaties onderling. Dit komt deels doordat de randomvariabele getrokken uit de LAS distributie soms zeer grote waarden kan aannemen. Het gevolg is dat de PDF van ∆x(t) dan een te vette staart krijgt enkel en alleen omdat een zeer onwaarschijnlijk event toevallig binnen onze beperkte simulatietijd valt. Een oplossing die hier gebruikt is, is een maximum op te leggen aan de waarde van de LAS-randomvariabele. Een betere oplossing is natuurlijk veel langer simuleren. Een andere vraag is dan ook of de vorm van de CDF van ∆x(t) in figuur 5.18 blijft behouden voor langere simulaties. Een andere opmerking is dat het in de eerste plaats eigenlijk de bedoeling was om aan de hand van de entropie (en de volatiliteit) een gelijkenis aan te tonen tussen de dynamica van


63

NEMD simulatie en die van de returns van de DAX index. Een goede vergelijking tussen de twee is echter slechts mogelijk als we een goede mapping hebben tussen markttijd een simulatietijd. In sectie 3.2.3.6 hebben we met behulp van de snelheidsautocorrelatiefunctie een ruwe schatting gemaakt, namelijk dat 1 dag marktijd overeenkomt met 300 simulatiestappen. Dit was echter voor een evenwichts MD simulatie. In figuur 5.19 is de snelheidscorrelatiefunctie van een NEMD simulatie gegeven. Uit deze figuur volgt dan de correspondentie 1 marktdag = 150 simulatiestappen. Voor ander NEMD simulaties kan dit nog sterk veranderen. De correlatie kan al op nul vallen na 30 stappen of pas na 500 stappen. Aangezien de mapping van markttijd en simulatietijd hier met de natte vinger is gedaan en alles van 50 tot 500 tijdstappen per marktdag een even juiste keuze is, kunnen we de entropie van de DAX index niet vergelijken met die van de NEMD simulatie. De entropie is bij de DAX index berekent over een periode van 200 dagen. De periode om de entropie in de simulatie te berekenen kan dan lopen van T = 200 × 50 = 10.000 tot T = 200 × 500 = 10.000. Een simulatie die amper

100.000 tijdstappen bevat is dus te kort. De keuze die we hier maken heeft een grote invloed op het verloop van de entropie, te zien in figuur 5.20. Als men tijdsafhankelijke grootheden wil vergelijken is het essentieel om de correcte tijdsmapping te hebben. Omdat we hier niet zeker waren van de juistheid van de mapping is er voor gekozen de vergelijking tussen NEMD simulatie en marktindex te maken aan de hand van de PDF’s en de CDF’s van ∆x(t). Figuur 5.19: De snelheidscorrelatiefunctie van een NEMD simulatie met temperatuursverloop weergegeven in 5.16.


64

Figuur 5.20: De entropie van de NEMD simulatie met temperatuursverloop weergegeven in 5.16, berekend met verschillende waarden voor T . De keuze van T heeft een grote invloed op het karakter van de entropie.

Hoofdstuk 6

Conclusie Het doel van deze masterproef was om de dynamica van de markt te beschrijven aan de hand van een moleculaire dynamica simulatie die op een of ander manier uit evenwicht wordt gebracht. In hoofdstuk 2 werden de returns, de volatiliteit, de informatiëentropie en de PDF van de returns van de DAX index in de periode 1991-2011 vergeleken met die van eenvoudig Geometrische Brownse Beweging model. De bevinding was dat het GBB model de complexe markten niet op een correcte manier beschrijft. Het grootste gebrek in de GBB simulatie was dat zeer grote en lage returns veel minder vaak voorkomen dan in een echte markt. Er werd aangetoond in figuur 2.10 dat de staarten van de PDF van de returns van de DAX index een machtswet volgen. Een dergelijke machtswet moet dus ook zeker in ons model voorkomen. Gewapend met de kennis van de dynamica van de DAX index zochten we een manier om die dynamica na te bootsen in ons model. Om het moleculair systeem uit evenwicht te brengen werd er geopteerd dit te doen door de temperatuur te veranderen (een andere manier is bv. de dichtheid te veranderen (zie bv. [20] en [21]). Na veel ’trail’ en ’error’ werd er een manier gevonden om het systeem zo uit evenwicht te brengen dat we in positieve staart van de PDF van ∆x(t) ongeveer een machtswet zagen (zie figuur 5.18). Alhoewel dit een bevredigend resultaat is, moeten we er wel een groot vraagteken bij zetten. Eerst en vooral moet er nog onderzocht worden in welke mate dit resultaat niet een toevalstreffer was. Om meer zekerheid te verkrijgen is het nodig dat de simulatie veel langer moet lopen om statistisch relevante resultaten te bekomen. Pas als we dan nog steeds een PDF uitkomen waarvan de positieve staart een machtswet vormt kunnen we stellen dat ons model geslaagd is.

65

Bibliografie [1] http://finance.yahoo.com. [2] L. Bachelier. Théory de la spéculation. Annales scientifiques de l’E.N.S 3e série, 17:21– 86, 1900. [3] Myron Scholes Fischer Black. The pricing of options and corporate liabilities. The Journal of Political Economy, 81(3):637–656, 1973. ¨ [4] A. Einstein. Uber die von der molekularkinetischen theorie der wärme geforderte bewegung von in ruhenden fl¨ ussigkeiten suspendierten teilchen.

Annalen der Physik,

322(8):549–560, 1905. [5] J. W. Lindeberg. Eine neue herleitung des exponentialgesetzes in der wahrscheinlichkeitsrechnung. Mathematische Zeitschrift, 15(1):211–225, 1922. [6] A. N. Kolmogorov V. V. Gnedenko. Limit Distributions of Sums of Independent Random Variables. 1954. [7] C.E. Shannon. A mathematical theory of communication. The Bell System Technical Journal, 27:623–656, 1948. [8] J. Voit. Statistical Mechanics of Financial Markets. Springer, 3 edition, 2005. [9] Jean-Philippe Bouchaud. Economics need a scientific revolution. Nature, vol 455, p. 1181 (2008), 2008. [10] H. Eugene Stanley Rosario N. Mantegna. Introduction to econophysics. Cambridge University Press, 2005. [11] Neil Johnson, Guannan Zhao, Eric Hunsader, Jing Meng, Amith Ravindar, Spencer Carran, and Brian Tivnan. Financial black swans driven by ultrafast machine ecology, 2012.

66

67

Bibliografie

[12] M. E. J. Newman. Pareto laws, pareto distributions and zipf’s law. Contemporary Physics, 46:323–351, 2005. [13] Nicholas R. Moloney Kim Christensen. Complexity and Criticality. Imperial College Press, 2005. [14] M. Waroquier. Kwantummechanica, februari 2009. Cursus nota’s. [15] J. Ryckebusch S. Standaerd and L. De Cruz. Creating conditions of anomalous selfdiffusion in a liquid with molecular dynamics. Journal of Statistical Mechanics, page P04004, 2010. [16] G. Franzese. Differences between discontinuous and continuous soft-core attractive potentials: The appearance of density anomaly. Journal of molecular liquids, 136:267–373, 2007. [17] R. Morent. Statistische fysica. Slides Academiejaar 2008-2009. [18] J. Ryckebush. Statistische fysica i, 2009. [19] H. Gould and J. Tobochnik. Thermal and Statistical Physics: With Computer Applications. Princeton University Press, 2007. [20] L. Baudewyn.

Simulatie van de volatiliteit in financiële markten met behulp

van niet-evenwichts moleculaire dynamica.

http://inwpent5.ugent.be/papers/

LiesbethBaudewyn.pdf. masterproef Academiejaar 2009-2010. [21] P. Van Nuffel.

Simulatie van informatie-entropie in financiëmarkten op basis

van niet-evenwichts moleculaire dynamica.

http://inwpent5.ugent.be/papers/

PieterVannuffel.pdf. masterproef Academiejaar 2010-2011.

Lijst van figuren 2.1

De logaritmische returns (∆t = 1 dag) van de BEL20 index in de periode van juli 2005 tot begin 2011. Om de returns te berekenen is de slotkoers van elke beursdag gebruikt [1]. De returns vertonen op bepaalde momenten uitgesproken pieken die overeenkomen met grote winsten en verliezen (zie bv de crash in september 2008). Deze pieken, en het feit dat ze geclusterd zijn, worden ook voor bijna elke andere index teruggevonden. . . . . . . . . . . . .

2.2

4

Volatiliteit van de BEL20 index in de periode van juli 2005 tot begin 2011. De volatiliteit kan gezien worden als de nervositeit of de ’temperatuur’ van de markten. Dit is duidelijk te zien in het tweede deel van het jaar 2008 door de financiële crisis in de tweede helft van 2008. . . . . . . . . . . . . . . . . . . .

2.3

5

In deze figuur is de informatieëntropie (Hdisc − ln ∆) weergegeven in functie van het aantal bins. De x-as is in logaritmische schaal. De entropie stijgt in

2.4

functie van het aantal bins (∆ → 0) . . . . . . . . . . . . . . . . . . . . . . .

10

De figuur toont een GBB simulatie met verschillende parameters (∆t = 0, 01 en 10000 tijdstappen). De parameter µ bepaalt hoe snel de prijs stijgt, σ bepaalt daarintegen hoe sterk de prijs fluctueert rond de exponentiële groei. . . . . .

2.5

12

Een GBB simulatie met parameters bepaalt uit het gemiddelde en de standaardafwijking van de te simuleren SP500 index in de periode 1991-2012(µ = 0, 0001565, σ = 0, 0052). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2.6

13

Een GBB simulatie met parameters bepaalt uit het gemiddelde en de standaardafwijking van de te simuleren NIKKEI index in de periode 1991-2012

2.7

(µ = −0, 000076213, σ = 0, 0067). . . . . . . . . . . . . . . . . . . . . . . . . .

14

Het bepalen van de parameters voor de GBB simulatie kan gebeuren door: 1. Een Gaussische functie te fitten (groen) aan de PDF van de returns (rood) of 2. Door met het gemiddelde en de standaarafwijking van de returns te berekenen (blauw). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

68

15

Lijst van figuren

2.8

69

Een vergelijking van de genormaliseerde returns van de DAX index en een GBB. De DAX index (groen) samen met een GBB simulatie met de random variabele getrokken uit een gaussische distributie (rood) en uit een Lévy alpha stable distributie (µ = 0, 5, a = 1, β = 0, m = 0) (blauw) . . . . . . . . . . . .

2.9

16

De probabiliteitsdistributiefuncties (PDF’s) voor de slotkoersen van de DAX index (groen) in de periode 1991-2011, een GBB met gaussische randomvariabele (rood) en een GBB met LAS random variabele (µ = 0, 5, a = 1, β = 0, m = 0) (blauw) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

17

2.10 De cumulatieve distributiefuncties (CDF) voor de beide helften ([−10, 0[ en [0, 10]) van de PDF (figuur 2.9) van de DAX index. De staarten van de CDF’s zijn gefit aan een machtswet y = Cxα . De exponent van de fit aan de positieve returns is −3, 73, die van de fit aan de negatieve −5, 94. . . . . . . . . . . . .

18

de DAX (groen) index. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

19

2.11 Een vergelijking van de genormaliseerde volatiliteit van de GGB (rood) en van 2.12 Vergelijking genormaliseerde Informatieëntropie GBB (rood) en DAX index (groen). De PDF’s van de genormaliseerde returns werden berekend in het interval [−10, 10] met bingrootte 0, 02. . . . . . . . . . . . . . . . . . . . . . .

20

2.13 Entropie (groen)en volatiliteit (rood) van de DAX index in de periode 19912011. De dynamica van de entropie is zeer gelijkaardig aan die van de volatiliteit. 21 2.14 De genormaliseerde informatieëntropie van de DAX index (groen) berekend met 2 bins vergeleken met de entropie van een Gaussische distributie. Zoals voorspeld is de entropie van de Gaussische distributie steeds maximaal. . . .

22

2.15 Autocorrelatiefunctie van de DAX index in de periode 1991-2011. De returns (blauw) zijn niet gecorreleerd aangezien de correlatiefunctie direct naar nul valt. De volatiliteit (rood) daarintegen vertoont wel een lange tijds correlatie die maar zeer traag daalt. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

23

2.16 Autocorrelatiefunctie van de GBB simulatie. De returns (blauw) zijn niet gecorreleerd. Deze eigenschap zit impliciet in het stochastisch proces. Het Wiener proces is immers een Markov proces. De volatiliteit vertoont een grotere correlatie maar valt redelijk snel lineair naar nul. . . . . . . . . . . . . . . . . 3.1

Een vergelijking van de Lennard-Jones potentiaal en de Softcore potentiaal met parameters uit tabel 3.1. De afstand r staat uitgedrukt in systeemeenheden. .

3.2

24

28

Het verloop van de temperatuur tijdens een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. Zoals verwacht blijft de temperatuur ongeveer dezelfde in dit gesloten systeem. . . .

31

Lijst van figuren

3.3

De energieën van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. . . . . . . . . . . . . . . . . . . . . . .

3.4

70

32

Het verloop van de gemiddelde verandering van de deeltjes van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.5

33

Probabiliteitsdistributiefunctie ∆x van de deeltjes van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. Zoals verwacht is dit in evenwicht een gaussische distributie. . . . . . . . . . .

3.6

34

Entropie in functie van de tijd van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. Zoals verwacht uit de thermodynamica schommelt de entropie een weinig rond een constante waarde voor een systeem in evenwicht. . . . . . . . . . . . . . . . . . . . . . .

3.7

Volatiliteit in functie van de tijd van een MD simulatie met 256 deeltjes en een initiële temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. . . . . . . . . . .

3.8

35 36

Snelheidsautocorrelatiefunctie van een MD simulatie met 256 deeltjes en een initiëele temperatuur T ∗ = 1, 1 en initiële druk van ρ∗ = 0, 5. . . . . . . . . .

37

4.1

Het kanonisch ensemble (figuur uit [17]) . . . . . . . . . . . . . . . . . . . . .

40

5.1

Het temperatuursverloop van een NEMD simulatie waarbij op tijdstappen 20.000, 40.000, 60.000 en 80.000 de temperatuur van het systeem plots wordt verhoogd of verlaagd. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.2

47

Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.1. De twee pieken stellen de ∆x’s voor in de twee evenwichtssituaties op temperatuur T ∗ = 0, 7 en T ∗ = 0, 95. Ten gevolge van de temperatuurschokken nemen we ook extremere waarden waar voor bepaalde ∆x. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.3

48

De entropie en de volatiliteit in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.1. De stijging en daling in entropie en volatiliteit zijn het gevolg van de plotse temperatuursveranderingen. . . . . .

5.4

49

Het temperatuursverloop van een NEMD simulatie waarbij tussen tijdstappen 10.000 en 20.000 de temperatuur van het systeem geleidelijk aan wordt verdubbeld. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

50

71

Lijst van figuren

5.5

Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.4. De twee pieken stellen de ∆x’s voor in de twee evenwichtssituaties op temperatuur T ∗ = 0, 7 en T ∗ = 1, 3. Ten gevolge van de temperatuursgradiënt worden er waarden van ∆x’s waargenomen die zich tussen de twee evenwichtspieken bevinden. . . . . . . . . . . . . . . . . . . . .

5.6

51

De entropie en de volatiliteit in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.4. De kleine stijging en daling in entropie en volatiliteit zijn het gevolg van de geleidelijke temperatuursveranderingen. .

5.7

52

Het temperatuursverloop van een NEMD simulatie waarbij na elke 2.000 tijdstappen de temperatuur gedurende 200 tijdstappen wordt verhoogd of verlaagd met een grootte bepaald door een gaussische randomvariabele.

5.8

. . . . . . . .

53

Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.7. De asymmetrie in het histogram is een gevolg van de implementatie van de temperatuurschok. . . . . . . . . . . . . . . . . .

5.9

54

De cumulatieve distributiefunctie van de positieve staart van het histogram in figuur 5.8 vergeleken met de CDF van een GBB simulatie en de CDF van de DAX index. Alhoewel grotere ∆x meer waarschijnlijk zijn, vertoond de positieve staart geen verval volgens een machtswet, maar vervalt eerder exponentieel. 55

5.10 De cumulatieve distributiefunctie van de negatieve staart van het histogram in figuur 5.8 vergeleken met de CDF van een GBB simulatie en de CDF van de DAX index. De negatieve staart vertoont geen van de karakteristieken die we zoeken. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

55

5.11 De entropie en de volatiliteit in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.1. De stijging en daling in entropie en volatiliteit zijn het gevolg van de plotse temperatuursveranderingen. . . . . .

56

5.12 Het temperatuursverloop van een NEMD simulatie waarbij na elke 2.000 tijdstappen de temperatuur gedurende 200 tijdstappen wordt verhoogd of verlaagd met een grootte bepaald door een randomvariabele getrokken uit een LAS distributie met parameters a = 0, 01, β = 0, µ = 0, 5 en m = 0. . . . . . . . . . .

57

5.13 Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.7. De asymmetrie in het histogram is een gevolg van de implementatie van de temperatuurschok. . . . . . . . . . . . . . . . . .

58

5.14 De cumulatieve distributiefunctie van de positieve staart van het histogram in figuur 5.8 vergeleken met de CDF van een GBB simulatie en de CDF van de DAX index. Alhoewel er een deel van de CDF vervalt als een machtswetswet, is dit niet het resultaat dat we zoeken. . . . . . . . . . . . . . . . . . . . . . .

59

Lijst van figuren

72

5.15 De cumulatieve distributiefunctie van de negatieve staart van het histogram in figuur 5.8 vergeleken met de CDF van een GBB simulatie en de CDF van de DAX index. Alhoewel de CDF vervalt als een machtswetswet, zijn de grote ∆x’s veel minder waarschijnlijk dan bij de CDF van de returns van de DAX index. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

59

5.16 Het temperatuursverloop van een NEMD simulatie waarbij na elke 1.000 tijdstappen de temperatuur gedurende 200 tijdstappen wordt verhoogd of verlaagd met een grootte bepaald door een randomvariabele voor 30% getrokken uit een LAS distributie met parameters a = 0, 01, β = 0, µ = 0, 5 en m = 0 en voor 70% getrokken uit een gaussische distributie met µ = 0 en σ = 1. . . . . . . .

60

5.17 Het histogram van ∆x(t) in een NEMD simulatie waarvan het temperatuursverloop gegeven is in figuur 5.16. De asymmetrie in het histogram is een gevolg van de implementatie van de temperatuurschok. . . . . . . . . . . . . . . . . .

61

5.18 De cumulatieve distributiefunctie van de positieve staart van het histogram in figuur 5.8 vergeleken met de CDF van de DAX index. De CDF vertoont grote gelijkenissen met die van de DAX index. . . . . . . . . . . . . . . . . . . . . .

62

5.19 De snelheidscorrelatiefunctie van een NEMD simulatie met temperatuursverloop weergegeven in 5.16. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

63

5.20 De entropie van de NEMD simulatie met temperatuursverloop weergegeven in 5.16, berekend met verschillende waarden voor T . De keuze van T heeft een grote invloed op het karakter van de entropie. . . . . . . . . . . . . . . . . . .

64

Lijst van tabellen 3.1

Parameters (in systeemeenheden) van gefitte SC potentiaal [15] . . . . . . . .

73

27

Niet-evenwichts moleculaire dynamica als een model voor financiële markten

Recommend Documents