Importance Sampling for Least Square Monte Carlo [duplicate]

Question

I am currently trying to implement and model an Importance Sampling estimator for Longstaff and Schwartz algorithm for pricing American put options. It is used such that more paths are in-the-money such that less simulation is required and also the variance is reduced.

I am following the steps by Moreni and simulate $n$-paths of a drifted geometric brownian motion given by:

$dS_t = S_t\big[(r+\theta \sigma)dt+\sigma dW_t\big]$.

Then moreni defines the likelihood ratio $L_\tau^\theta=exp\{-\theta W_t+0.5 \theta^2 t \}$.

By the stopping time theorem, we want to find for each path $n$ $\sup_\tau \mathbb{E}[f^{(n)}(\tau,S_\tau)]$ the corresponding Importance sampling estimator $\sup_\tau \mathbb{E}[L_\tau^\theta f^{(n)}(\tau,S_\tau)]$, where $f$ is the payoff function.

When we have found the stopping time for each path $n$ we can take the average and get the desired put price at time 0.

My question is that is the wiener process in the likelihood ratio the random number at $\textit{optimal}$ time $\tau$ for path some path $n_i$? And at the same time should $t$ be the optimal time $\tau$ for path $n_i$?

Since in R the standard L&S algorithm is existent in the package library(LSMonteCarlo), i modified the function AmerPutLSM, such that it simulate with drift $\theta$ and i set the interest rate to 0. However, I seem to get very high prices when I try to simulate compared to the standard one. So I am not sure if I really have understood the math correctly. I use negative values for $\theta$ and it should be around -0.5 and -1 according to Moreni. But the prices I get are really deviating.

Here is the code, that I modified:

firstValueRow <- function (x)
{
  cumSumMat <- matrix(NA, nrow = dim(x)[1], ncol = dim(x)[2])
  for (i in 1:(dim(x)[1])) {
    cumSumMat[i, ] <- cumsum(x[i, ])
  }
  cumSumMat2 <- cbind(matrix(0, nrow = dim(x)[1], ncol = 1), cumSumMat[, -(dim(x)[2])])
  ResultMat <- matrix(NA, nrow = dim(x)[1], ncol = dim(x)[2])
  for (i in 1:dim(x)[2]) {
    ResultMat[, i] <- ifelse(cumSumMat2[, i] > 0, 0, x[,i])
  }
  return(ResultMat)
}

Spot = 36
sigma = 0.8
theta = -1
mu = 0
n = 1000
m = 50
Strike = 40
r = 0
dr = 0
mT = 1

dt <- mT/m
GBM <- matrix(NA, nrow = n, ncol = m)
Zlist <- matrix(NA, nrow = n, ncol = m)
for (i in 1:n)
{
  Z <- rnorm(m, mean = mu, sd = 1)
  GBM[i, ] <- Spot * exp(cumsum(((r+theta*sigma-0.5*sigma*sigma)*dt)+(sigma*(sqrt(dt))*Z)))
  Zlist[i, ] <- Z
}
X <- ifelse(GBM < Strike, GBM, NA) #stock rates only in the money
### payoff importance sampling
Ltheta <- exp(-theta*Zlist-0.5*theta*theta*dt)
CFL <- matrix(Ltheta*pmax(0, Strike - GBM), nrow = n, ncol = m) #cashflows

Xsh <- X[, -m]
X2sh <- Xsh * Xsh

Y1 <- CFL * exp(-1 * r * dt)
Y2 <- cbind((matrix(NA, nrow = n, ncol = m - 1)), Y1[, m]) #value of derivate at time t+1

CV <- matrix(NA, nrow = n, ncol = m - 1) #continuation value
for (i in (m - 1):1)
{
  reg1 <- lm(Y2[, i + 1] ~ Xsh[, i] + X2sh[, i])
  CV[, i] <- (matrix(reg1$coefficients)[1, 1]) + ((matrix(reg1$coefficients)[2,1]) * Xsh[, i]) + ((matrix(reg1$coefficients)[3,1]) * X2sh[, i])
  CV[, i] <- (ifelse(is.na(CV[, i]), 0, CV[, i]))
  Y2[, i] <- ifelse(CFL[, i] > CV[, i], Y1[, i], Y2[, i + 1] * exp(-1 * r * dt))
}
CV <- ifelse(is.na(CV), 0, CV)
CVp <- cbind(CV, (matrix(0, nrow = n, ncol = 1)))
POF <- ifelse(CVp > CFL, 0, CFL)
FPOF <- firstValueRow(POF)
dFPOF <- matrix(NA, nrow = n, ncol = m)
for (i in 1:m)
{
  dFPOF[, i] <- FPOF[, i] * exp(-1 * dt * r * i)
}

PRICE <- mean(rowSums(dFPOF))
PRICE

Should the cashflows be multiplied by the likelihood ratio immediately?

Appreciate for help. Thanks

Answer 1

I just made things clearer hoping it would help.

Let define $\mathbb{Q}_\theta$ as $$\frac{d\mathbb{Q}_\theta}{d\mathbb{P}}|_{\mathcal{F}_t}=\exp(\theta W_t -\frac{1}{2}\theta^2 t)=Z^\theta_t$$

By girsanov, if $W$ is a brownian motion under $\mathbb{P}$, then $W^\theta_t=W_t-\theta t$ is a brownian motion under $\mathbb{Q}^\theta$

$$\begin{split} \mathbb{E}[f(\tau, W_\tau)]=& \mathbb{E}^{\mathbb{Q}_\theta}[(Z_\tau^\theta)^{-1} f(\tau,W_\tau)] \\ = & \mathbb{E}^{\mathbb{Q}_\theta}[e^{-\theta W^\theta_\tau-\frac{1}{2}\theta^2\tau}f(\tau,W^\theta_\tau+\theta\tau)] \end{split}$$

and if $\tau$ was defined as :

$$\tau = \inf\{t\geq 0 : W_t \in A_t\}$$ then $$\tau = \inf\{t\geq 0 : W^\theta_t+\theta t \in A_t\}$$

defining $\tau^\theta$ as :

$$\tau^\theta= \inf\{ t\geq 0 : W_t+\theta t \in A_t\}$$

you finally get (using your $L$ notation) and the equality in law

$$\mathbb{E}[f(\tau, W_\tau)]=\mathbb{E}[L_{\tau^\theta}^\theta f(\tau^\theta,W_{\tau^\theta}+\theta \tau^\theta)]$$