R code for simulating stochastic asset price path

Question

Consider the following model for the evolution of an asset's price:

This what I have done (in R). I could not find a function that randomly outputs +1 or -1, so I decided to adapt the inbuilt rbinom function.

## This code is in R

rm(list = ls())

library(dplyr)
library(dint)
library(magrittr)
library(stats)

path  = 
  function(T, mu, sigma, p, x0) {
    x    = rep(NA, T)
    x[1] = x0
    
    for(i in 2:T){
    z    = if_else(rbinom(1,1,p) == 0, -1, 1)
    x[i] = x[i-1] * exp(mu + sigma*z)
    }
    return(x)
  }

## Just some testing
x_sim = path(T = 4, mu = 0, sigma = 0.01, p = 0.5, x0 = 100)

## Actual answer
Np = 10000
mc = matrix(nrow = 17, ncol = Np)
for(j in 1:Np){
  mc[,j] = path(T = 17, mu = 0, sigma = 0.01, p = 0.5, x0 = 100)
}
test     = mc[2:nrow(mc), ] >= 100
sum_test = colSums(test)
comp     = sum(sum_test >= 1)/length(sum_test)
prob     = 1 - comp

Does this make sense? Any help/tips/advice would be much appreciated. Thanks!

Answer 1

Staying close to your code, I came up with this. Intuitively, if you think about it, the probability should be rather low due to the parameters and I get a probability of about 6.7% which is roughly what I get if I run your code with the parameters from the assignment.

simpath <- function(t, mu, sigma, p, x0, seed){
  # set seed
  if(!missing(seed)){
    set.seed(seed)  
  }
  # set up matrix for storing the results
  res <- matrix(c(1:t, rep(NA, t*2)), ncol = 3)
  colnames(res) <- c('t', 'z_t', 'x_t')
  res[, 'z_t'] <- sample(c(1, -1), size = t, prob = c(p, 1-p), replace = TRUE)
  res[1, 3] <- x0
  for(i in 2:t){
    res[i, 3] <- res[i-1, 3] * exp(mu+sigma*res[i, 2])
  }
  return(res)
}

x_sim <- simpath(t = 4, mu = 0, sigma = 0.01, p = 0.5, x0 = 100, seed = 123)
x_sim2 <- simpath(t = 36, mu = 0, sigma = 0.03, p = 0.5, x0 = 100, seed = 123)

## Actual answer
Np <- 100000
mc <- matrix(nrow = 36, ncol = Np)
for (j in 1:Np){
  mc[, j] <- simpath(t = 36, mu = 0, sigma = 0.03, p = 0.5, x0 = 100)[, 3]
}
test <- mc > 100
sum_test <- colSums(test)
comp     = sum(sum_test == 0)/length(sum_test)
prob     = comp
> prob
[1] 0.06759