R/Merton.R

In gaspardcc/jumpyMerton:

#' Merton class
#'
#' Inheritance from Underlying class
#'
#' @author elh

Merton <- setClass(
  # Name of the class
  "Merton",

  slots = c(
    jumpAvg = "numeric",
    jumpVar = "numeric",
    jumpIntensity = "numeric"
  ),

  prototype = list(
    jumpAvg = 0.05,
    jumpVar = 0.2,
    jumpIntensity = 0.5
  ),

  # Inheritance
  contains = "BlackScholes"
)

setMethod(
  f = "asset",
  signature = "Merton",
  definition = function(this, maturity, timestep)
  {
    today <- tail(this@grid, n=1)

    if(length(today)==0)
      today = 0
    # number of jump
    njump <- rpois(1, this@jumpIntensity * maturity)

    # Place of those jumps into the timeline
    jump <- runif(njump, today, today+maturity)

    # Set the grid with jumps
    this@grid <- sort(c(jump, seq(from = today, to=today+maturity, by=timestep)))

    # Specify the jump
    Y <- rnorm(njump, this@jumpAvg, this@jumpVar)

    # Brownian increment
    G <- rnorm(maturity/timestep + njump + 1)

    # Logarithm of the the underlying path
    X <- vector("numeric", length(this@grid))
    if(length(this@path)==0)
      X[1] <- 0
    else
      X[1] <- log(tail(this@path, n=1)/this@spot)

    for(i in 2:length(this@grid))
    {
      delta_t = this@grid[i] - this@grid[i-1]
      X[i] <- X[i-1] + this@drift * delta_t + this@volatility * sqrt(delta_t) * G[i]

      if(this@grid[i] %in% jump) X[i] <- X[i] + Y[which(this@grid[i] == jump)]
    }
    this@path <- this@spot * exp(X)
    return(this)
  }
)

setMethod(
  f = "phi",
  signature = "Merton",
  definition = function(this, v, maturity)
  {
    return(exp(maturity * (-this@volatility^2/2 + complex(imaginary = 1)*v + this@jumpIntensity*(exp(-(this@jumpAvg*v)^2/2+ complex(imaginary = 1)*this@drift*v)-1))))
  }
)