Opakovanie: Wienerov proces

Wienerov proces je súbor náhodných premenných s parametrom \(t\) (čas), pre ktorý platí:

Úloha 1: Rozdelenie \(w_1 + w_2\)

Aké rozdelenie má náhodná premenná \(w_1 + w_2\)?

Úloha 2: Simulácia Wienerovho procesu

Naprogramujeme funkciu, ktorá vráti hodnoty Wienerovho procesu v \(n\) dátových bodoch vzdialených \(\Delta t\):

Wiener <- function(n, delta_t) {
  w <- numeric(n)
  
  w[1] <- 0
  for (i in 2:n) {
    w[i] <- w[i-1] + rnorm(1, mean = 0, sd = sqrt(delta_t))
  }
  
  return(w)
}

Rozmyslite si, ako táto funkcia funguje.

Vykreslíme viacero trajektórií. Vašou úlohou je opraviť hodnotu ylim tak aby koncový bod trajektórie bol na \(95\%\) v grafe.

n <- 201
delta_t <- 0.005

t <- seq(0, (n - 1)*delta_t, by = delta_t)
ylim <- c(-1, 1)  # FIXME

w <- Wiener(n, delta_t)
plot(t, w, type = "l", lwd = 2, col = rgb(0, 0, 0, 0.9), ylim = ylim)

for (i in 1:50) {
  w <- Wiener(n, delta_t)
  lines(t, w, col = rgb(0, 0, 0, 0.2))
}

Bonus: Celá trajektória (1b)

Upravte premennú ylim tak, aby celá trajektória bola na \(95\%\) v grafe. Návod: pozrite si proces \(M(t) = \sup \limits_{s \leq t} w(t)\), známy ako running maximum.

Bonus: Wienerov proces bez cyklu (1b)

Naprogramujte funkciu, ktorá robí to isté, ale neobsahuje for-cyklus. Porovnajte časy potrebné na zbehnutie týchto dvoch metód.

tic <- Sys.time()
w <- Wiener(1000000, 0.01)
toc <- Sys.time()

print(round((toc-tic), 5))
## Time difference of 1.82995 secs
WienerLoopless <- function(n, dt) {
  # TODO
}

tic <- Sys.time()
w <- WienerLoopless(1000000, 0.01)
toc <- Sys.time()

print(round((toc-tic), 5))
## Time difference of 0.00199 secs

Bonus: Zlomkový Brownov pohyb (2b)

Ak namiesto nezávislých prírastkov budeme uvažovať proces, ktorý spĺňa: \[ Cov(x(t), x(s)) = \frac12 \left(t^{2H} + s^{2H} - |t-s|^{2H}\right), \] pre parameter \(H \in (0, 1)\), dostaneme tzv. zlomkový Brownov pohyb. Parameter \(H\) sa nazýva Hurstov index.

Napíšte funkciu, ktorá bude generovať trajektórie takéhoto procesu. Vykreslite 10 trajektórii postupne pre \(H = 0.1, 0.25, 0.5, 0.75, 0.9\). Slovne popíšte aký efekt má parameter \(H\) na priebeh procesu. použite n = 1001 časových bodov s rozostupmi dt = 0.01.

Návod: vytvorte kovariančnú maticu \(\Sigma\) a vypočítajte jej maticovú odmocninu \(\Gamma\), t.j. \(\Sigma = \Gamma^T \Gamma\). (Napríklad pomocou spektrálneho rozkladu matice \(\Sigma\).) Potom vytvorte vektor \(v\), ktorý obsahuje nezávisle generované hodnoty z normálneho rozdelenia so strednou hodnotou \(0\) a disperziou \(1\). Napokon hodnoty vektora \(\Gamma v\) predstavujú trajektóriu zlomkového Brownovho pohybu.

FractionalBrownianMotion <- function(n, dt) {
  # TODO
}

Úloha 3: \(w_1 + w_2\) numericky

Vygenerujte hodnoty náhodnej premennej \(w_1+w_2\) a spravte histogram jej hodnôt. Všimnime si, že nám stačí generovať Wienerov proces v časoch \(1\) a \(2\), jemnejšie delenie nie je potrebné. Porovnajte pravdepodobnostné rozdelenie tejto náhodnej premennej so získaným histogramom.

# TODO

Úloha 4: Brownov pohyb (vychýlený Wienerov proces)

Nech \(w_t\) je Wienerov proces. Potom proces \[ B_t := \mu t + \sigma w_t, \] kde \(\mu, \sigma\) sú konštanty, nazývame Brownovým pohybom.

Poznámka: Terminológia tu nie je jednoznačná. Niektorí nazývajú Brownovým pohybom Wienerov proces (bez driftu, s jednotkovou volatilitou). Iní zas tvrdia, že Brownov pohyb je fyzikálny fenomén a nie stochastický proces. Pri práci s literatúrou si vždy pozrite značenie.

Naprogramujte funkciu BrownianMotion(n, delta_t, mu, sigma) tak, aby vygenerovala trajektóriu Brownovho pohybu.

BrownianMotion <- function(n, delta_t, mu, sigma) {
  # TODO
}

Podobne ako minule, vykreslíme viac trajektórií do toho istého obrázku. Nastavte vertikálnu os tak, aby koncový bod trajektórie bol na \(95\%\) vnútri grafu.

n <- 201
delta_t <- 0.005

t <- seq(0, (n - 1)*delta_t, by = delta_t)
ylim <- c(-1, 1)  # FIXME

b <- BrownianMotion(n, delta_t, mu, sigma)
plot(t, b, type = "l", lwd = 2, col = rgb(0, 0, 0, 0.9), ylim = ylim)

for (i in 1:50) {
  b <- BrownianMotion(n, delta_t, mu, sigma)
  lines(t, b, col = rgb(0, 0, 0, 0.2))
}

Úloha 5: Geometrcký Brownov pohyb

Uvažujme náhodný proces \[ S_t = S_0 \exp(\mu t + \sigma w_t). \] Tento proces sa nazýva geometrický Brownov pohyb.

Naprogramujte funkciu, ktorá vráti hodnoty trajektórie geometrického Brownovho pohybu

GeometricBrownianMotion <- function(n, delta_t, mu, sigma, s0) {
  # TODO
}



n <- 201
delta_t <- 0.005

t <- seq(0, (n - 1)*delta_t, by = delta_t)
ylim <- c(-1, 1)  # FIXME

b <- BrownianMotion(n, delta_t, mu, sigma)
plot(t, b, type = "l", lwd = 2, col = rgb(0, 0, 0, 0.9), ylim = ylim)

for (i in 1:50) {
  b <- GeometricBrownianMotion(n, delta_t, mu, sigma)
  lines(t, b, col = rgb(0, 0, 0, 0.2))
}

Poznámka: Existujú dva štandardné tvary geometrického Borwnovho pohybu. Jeden z nich je \[ S_t = S_0 \exp(\mu t + \sigma w_t), \] ktorý je riešením stochastickej diferenciálnej rovnice\[ dS_t = \left(\mu + \frac{\sigma^2}{2}\right)S_t \, dt + \sigma S \, dt \]

Druhý tvar je \[ S_t = S_0 \exp\left( \left(\mu - \frac{\sigma^2}{2} \right)t + \sigma w_t \right) \] daný stochastickou diferenciálnou rovnicou\[ dS_t = \mu S \,dt + \sigma S \, dw_t. \]

Bonus: Nezáporný proces (1b)

Ukážte že proces daný stochastickou diferenciálnou rovnicou \[ dr_t = (2 - r_t)\ dt + r_t\ dw_t \] s počiatočnou podmienkou \(r_0 = 1\) nadobúda iba nezáporné hodnoty.

Opakovanie: Itôva lema

Itôva lema nám pomôže nájsť tvar stochastickej diferenciálnej rovnice pre transformovaný proces. Ak proces \(X\) daný SDR \(dX = \mu \, dt + \sigma \, dw_t\) transformujeme dvakrát diferencovateľnou funkciou \(g(t, x)\), potom pre proces \(Y := g(t, X)\) platí: \[ dY = \left(\frac{\partial g}{\partial t} + \mu \frac{\partial g}{\partial x} + \frac{\sigma^2}{2} \frac{\partial^2 g}{\partial x^2} \right)\,dt + \sigma \frac{\partial g}{\partial x}\, dw_t \]

Príklad

Majme stochastický proces daný stochastickou diferenciálnou rovnicou: \[ dX_t = 0,\!98\, dt + 0,\!5 \,dw_t \] Jedná sa o vychýlený Wienerov proces.

Nájdite stochastickú diferenciálnu rovnicu pre proces \(S := e^X\). Procesu \(S\) sa hovorí aj geometrický brownov pohyb.

Riešenie:

Funkcia \(g\), ktorou transformujeme proces má tvar \(g(t, x) = e^x\). Jej parciálne derivácie sú: \[ \begin{align*} \frac{\partial g}{\partial t} &= 0 \\ \frac{\partial g}{\partial x} &= e^x \\ \frac{\partial^2 g}{\partial x^2} &= e^x \end{align*} \]

Teda: \[ dS = \left(0 + 0,\!98e^X + \frac{(\frac15)^2}{2} e^X\right)\,dt + \frac15 e^X dw_t \] Namiesto \(e^X\) budeme písať \(S\) a dostaneme: \[ \boxed{dY = S \,dt + \frac15 S\, dw_t} \]

Iné riešenie:

Použijeme Taylorov polynóm: \[ g(t+dt, X+dX) - g(t, x) = \frac{\partial g}{\partial t} dt + \frac{\partial g}{\partial x}dX + \frac12 \frac{\partial^2 g}{\partial x^2}(dX)^2 + o(dt) + o(dx^2) \] \[ \begin{align*} dS &= 0\,dt + e^X dX + \frac12e^X(dX)^2 \\ &= e^X (0,\!98\, dt + 0,\!2 \,dw_t) + \frac12 e^X(0,\!98\, dt + 0,\!2 \,dw_t)^2 \\ &= S (0,\!98\, dt + 0,\!2 \,dw_t) + \frac12 S(0,\!98\, dt + 0,\!2 \,dw_t)^2 \\ &= S (0,\!98\, dt + 0,\!2 \,dw_t) + \frac12 S(\underbrace{0,\!9604\, dt^2}_{o(dt)} + \underbrace{0,\!392 \, dt \, dw}_{o(dt)} + 0,\!04 \,\underbrace{dw_t^2}_{dt}) \\ &= 0,\!98 S \, dt + 0,\!2 S\,dw_t + 0,\!02 S \,dt \\ &= \boxed{S \, dt + 0,\!2 S \, dw_t} \end{align*} \]

Úloha 6: Itôva lema

\[ dr_t = (2-r_t)\ dt + 3 \sqrt{r_t} \ dw_t \] S týmto stochastickým procesom sa ešte stretneme pri modelovaní úrokových mier.

Nájdite stochastickú diferenciálnu rovnicu pre proces \(f := r - r_0e^{-t} - 2\left(1-e^{-t}\right)\).

Vo výsledku sa nesmie objaviť hodnota procesu \(r\).

Riešenie (začiatok)

Funkcia, ktorou transformujeme má tvar

\[ \begin{align*} g(t, x) &= x - r_0 e^{-t} - 2\left( 1 - e^{-t} \right) \\ &= x + (2 - r_0)e^{-t} - 2 \end{align*} \]

Úloha 7: Itôva lema ešte raz

Majme stochastický proces daný stochastickou diferenciálnou rovnicou:

\[ dX_t = (3X_t - X_t^3)dt + dw_t \] Nájdite diferenciál pre proces \(Y := X^2\)

Generovanie trajektórií (Euler-Maruyamova metóda)

Namiesto infinitezimálneho časového kroku \(dt\) použijeme konečný časový krok \(\Delta t\). Namiesto infinitezimálneho prírastku procesu \(dX\) budeme uvažovať konečnú diferenciu \(X_{t + \Delta t} - X_t\) a namiesto infinitezimálneho prírastku Wienerovho procesu použijeme náhodnú premennú so strednou hodnotou \(0\) a disperziou \(\Delta t\).

Namiesto rovnice v tvare \[ dX_t = \mu(X_t, t)\, dt + \sigma(X_t, t) \, dw_t \] budeme počítať \[ X_{t+\Delta t} - X_t = \mu(X_t, t)\cdot \Delta t + \sigma(X_t, t)\cdot \mathcal N(0, \Delta t), \] resp. \[ X_{t+\Delta t} = X_t + \mu(X_t, t)\cdot \Delta t + \sigma(X_t, t)\cdot \mathcal N(0, \Delta t), \]

Poznámka: existujú aj iné numerické schémy na simuláciu náhodných procesov, napríklad Milsteinova metóda alebo metóda Leimkuhlera a Matthewsa.

Preštudujte si nasledujúci kód a vysvetlite prečo

EulerMaruyama <- function(Drift, Volatility, t_points, x0 = 0) {
  n <- length(t_points)
  dt <- diff(t_points)
  dw <- rnorm(n - 1, mean = 0, sd = sqrt(dt))
  
  x <- numeric(n)
  x[1] <- x0
  for (i in 1:(n-1)) {
    x[i+1] <- x[i] + 
              Drift(x[i], t_points[i]) * dt[i] + 
              Volatility(x[i], t_points[i]) * dw[i]
  }
  
  return(x)
}

Riešme pre \(dX = X\,dt + \frac15\,dw_t\) s počiatočnou podmienkou \(X_0 = 10\) na intervale \([0, 4]\) a s časovým krokom \(\Delta t = 1/365\).

x0 <- 10
t_max <- 4
delta_t <- 1/365

Drift <- function(x, t) {
  return(1 * x)
}

Volatility <- function(x, t) {
  return(0.2 * x)
}


t_points <- seq(0, t_max, by = delta_t)

x <- EulerMaruyama(Drift, Volatility, t_points = t_points, x0 = x0)
plot(t_points, x, "l", xlab = "t", ylab = "X(t)")

Bonus: Milsteinova metóda (1b)

Implementujte Milsteinovu metódu. Je to iná metóda na generovanie trajektórií náhodných procesov. Vykreslite trajektóriu vami vybraného procesu.

Milstein <- function(Drift, Volatility, t_points, x0 = 0) {
  # TODO
}

Bonus: Zlá simulácia (1b)

Nájdite príklad takého stochastického procesu, že Euler-Maruyamova metóda môže dať kvalitatívne zlý priebeh (napríklad trajektória, ktorá má byť nezáporná bude obsahovať hodnotu menšiu ako 0).

Bonus: Trajektórie pomocou Itôvej lemy (1b)

Simulujte priebeh stochastického procesu \(Y_t = X_t^2\), kde \(dX = (3X - X^3)\,dt + dw_t\). Zvoľte počiatočnú podmienku \(X_0 = 1\) a časový krok \(\Delta t \leq 0,\!01\). Urobte to dvoma spôsobmi:

  • Najprv použite Euler-Maruyamovu metódu na proces \(Y_t\).
  • Potom vygenerujte priebeh procesu \(X_t\), použijúc tie isté hodnoty \(dw_t\) ako v predošlom bode. Takto získanú trajektóriu transformujte funkciou \(g(t, x) = x^2\). Porovnajte získané trajektórie.
EulerMaruyama <- function(Drift, Volatility, t_points, x0 = 0, dw = NULL){
  n <- length(t_points)
  dt <- diff(t_points)
  
  if (is.null(dw)) {
    dw <- rnorm(n - 1, mean = 0, sd = sqrt(dt))
  }
  
  x <- numeric(n) # pre-alokácia
  x[1] <- x0
  for (i in 1:(n-1)) {
    x[i+1] <- x[i] + 
              Drift(x[i], t_points[i]) * dt[i] + 
              Volatility(x[i], t_points[i]) * dw[i]
  }
  
  return(x)
}
PrepareParameters <- function(t_max = 1000, delta_t = 0.01) {
  t_points <- seq(0, t_max, by = delta_t)
  
  n <- length(t_points)
  dt <- diff(t_points)
  dw <- rnorm(n - 1, sd = sqrt(dt))
  
  return(dw)
}

# Spôsob 1
# TODO

# Spôsob 2
# TODO

Úloha 8: Dvojstabilný proces

Simulujte proces \(X_t\) daný SDR \(dX_t = (3X_t - X_t^3)dt + dw_t\). Zvoľte časový interval \([0, 1000]\), počiatočnú podmienku \(X_0 = 1\) a časový krok \(\Delta t \leq 0,\!01\). Vykreslite viacero simulácii a slovne popíšte čo vidíte.

t_max <- 1000
delta_t <- 0.01
t_points <- seq(0, t_max, by = delta_t)

Drift <- function(x, t) {
  0 # FIXME
}

Volatility <- function(x, t) {
  1
}


for (i in 1:5) {
  x <- EulerMaruyama(Drift, Volatility, t_points = t_points, x0=1)
  plot(t_points, x, "l", xlab = "t", ylab = "X_t")
}

Bonus: stredná hodnota a disperzia (6b)

Túto úlohu môžete riešiť až do konca semestra.

  • Nájdite strednú hodnotu procesu daného stochastickou diferenciálnou rovnicou \[ dX_t = (3X_t - X_t^3)dt + dw_t \]

    v čase \(t\) za podmienky že \(X_0\) bolo rovné \(x_0\): \[ \mathbb E[X_t | X_0 = x_0] \]

  • Nájdite limitu tejto strednej hodnoty pre \(t \to \infty\).

  • Nájdite disperziu procesu v čase \(t\) a jej limitu pre \(t \to \infty\).

Poznámka: táto úloha je veľmi ťažká. Stačí ak úspešne vyriešite aspoň jej časť.

Úloha 9: Cena akcie

Predpokladajme, že cena akcie sa riadi geometrickým Brownovým pohybom: \[ dS = 0,\!3S\, dt + 0,\!25 S \, dw_t. \] Dnešná cena akcie je \(150\).

Bonus: Zlé odvodenie (0,5b)

Zdôvodnite prečo pre proces daný SDR \[ dX_t = (3X_t - X_t^3)\, dt + dw_t \] s počiatočnou podmienkou \(X_0 = x_0\) neplatí diferenciálna rovnica pre strednú hodnotu v tvare \[ \frac{df(t)}{dt} = 3f(t) - f^3(t) \] s počiatočnou podmienkou \(f(0) = x_0\).

Ďalšie úlohy

Viac úloh nájdete v archíve cvičení.

Odovzdanie bonusov

Bonusové úlohy odovzdajte do 6 dní mailom s predmetom FD24 – nick, kde nick je vami zvolená prezývka skupiny.

Veľká bonusová domáca úloha: stratégie na reálnych dátach (max. 3b)

Pozrite si historické dáta spoločnosti Apple Inc.

Vytvorte portfólio, ktoré:

Strike Bid Ask
135.00 48.05 48.75
140.00 43.10 43.80
145.00 38.20 38.85
150.00 33.25 33.90
155.00 28.30 28.95
160.00 23.50 23.90
165.00 18.70 19.00
170.00 14.00 14.25
175.00 9.70 9.85
180.00 6.5 6.15
185.00 3.30 3.40
190.00 1.59 1.64
195.00 0.72 0.74
200.00 0.33 0.34
205.00 0.17 0.19
210.00 0.09 0.13
215.00 0.06 0.09
220.00 0.03 0.05
225.00 0.01 0.05
230.00 0.00 0.04

Hodnotenie:

Odovzdať treba:

Vaše riešenie odovzdajte mailom s predmetom FD24 – bonus do 06.03.2024.

Domáca úloha (5b)

Úlohu vypracujte samostatne alebo v skupinkách najviac 3 ľudí.

Vaše riešenie odovzdajte mailom s predmetom FD24 – nick , kde nick je vami zvolená prezývka skupiny. Riešenie treba odovzdať do 06.03.2024.

Úloha 1: Itôva lema pre diskrétne prírastky (2b)

Podľa Itôvej lemy platí \(dw_t^2 = dt\). Ukážeme, že podobné tvrdenie platí aj pre diskrétny časový krok, aspoň pre strednú hodnotu.

Označme \(\Delta w\) zmenu Wienerovho procesu za časový interval dĺžky \(\Delta t\): \(\Delta w_t = w_{t + \Delta t} - w_t\).

  • Aké je rozdelenie náhodnej premennej \((\Delta w)^2\)? Akú ma strednú hodnotu a disperziu?
  • Ukážte že \(\mathbb E\left[(\Delta w)^2\right]\) a \(\Delta t\) sú asymptoticky zhodné, t.j. \[ \lim_{\Delta t \to 0^+} \frac{\mathbb E\left[(\Delta w)^2\right]}{\Delta t} = 1 \]

Úloha 2: Ornstein-Uhlenbeckov proces

Uvažujme proces daný stochastickou diferenciálnou rovnicou: \[ dr_t = \kappa (\theta - r_t)\, dt + \sigma\, dw_t \]

2a) Stredná hodnota procesu (1b)

Aplikujme strednú hodnotu na stochastickú diferenciálnu rovnicu: \[ \begin{align*} \mathbb E[dr_t] &= \mathbb E[\kappa (\theta - r_t)\, dt + \sigma\, dw_t] \\ d\mathbb E[r_t] &= \kappa (\theta - \mathbb E[r_t])\, dt + \sigma \mathbb E[dw_t] \end{align*} \] Stredná hodnota \(dw_t\) je rovná nule. Ak označíme \(f(t) := \mathbb E[r_t]\), dostávame diferenciálnu rovnicu pre strednú hodnotu: \[ df(t) = \kappa (\theta - f(t))\, dt \] Vyriešte túto diferenciálnu rovnicu s počiatočnou podmienkou \(f(0) = r_0\). Aká je limitná stredná hodnota tohto procesu, t.j. \[ \lim_{t\to\infty} \mathbb E[r_t] \]

2b) Parametre procesu (2b)

Priraďte hodnoty parametrov k trajektóriám na obrázku:

  • \(\kappa = 20, \theta = 1, \sigma = 1\)
  • \(\kappa = 3, \theta = 1, \sigma = 1\)
  • \(\kappa = 20, \theta = 3, \sigma = 5\)
  • \(\kappa = 20, \theta = 3, \sigma = 1\)

Návod: Pomôžte si simuláciou pomocou metódy Euler-Marayuma.