Euler–Maruyama method

OUprocess.R

library(ggplot2)
library(dplyr)
#browseURL("https://en.wikipedia.org/wiki/Euler–Maruyama_method")

mu <- function(y, t, mod){
  mod$THETA * (mod$MU - y)
}

sigma <- function(y, t, mod){
  mod$SIGMA
}

dW <- function(delta_t){
  rnorm(1, 0, sqrt(delta_t))
} 

run_simulation <- function(modparam){
  T_INIT = 3
  T_END = 7
  N = 1000  # Compute at 1000 grid points
  DT = (T_END - T_INIT) / N
  TS = seq(T_INIT, T_END + DT, by=DT)
  #assert TS.size == N + 1
  Y_INIT = 0
  ys = numeric(length(TS))
  ys[1] = Y_INIT
  for(i in 1:(length(TS)-1)){
    t = T_INIT + i * DT
    y = ys[i]
    ys[i+1] = y + mu(y, t, modparam) * DT + sigma(y, t, modparam) * dW(DT)
  }
  data.frame(time=TS, y=ys)
}

modparam <- list(
  THETA = 0.7,
  MU = 1.5,
  SIGMA = 0.06  
)

results <- lapply(1:5, function(i){
  set.seed(i)
  df = run_simulation(modparam)
  df$group = i
  df
})

result <- bind_rows(results)
ggplot(result, aes(x=time, y=y, group = group))+
  geom_line()+
  theme_bw(18)

ggsave("OUP.png")