<a href="https://colab.research.google.com/github/Jun-629/20MA573/blob/master/src/Hw11_IS.ipynb" target="_parent"></a>

Asset price under $\mathbb Q$ follows
$$S_{t} = S_{0} exp \{ \mu t + \sigma W_{t} \}$$
Coonsider Digital put with payoff 
$$h(S_{T}) = I(S_{T} < S_{0}e^{-b})$$
We want to find the forward price:
$$v = \mathbb E^{\mathbb Q}[h(S_{T})]$$
Parameters are given as 
$$r = .03, \sigma = .2, \mu = r - \frac{1}{2} \sigma^2 = .01, T = 1, b = .39$$


- Prove that the exact price is 0.02275

__Pf:__

As we have discussed in the video, we claim that

\begin{equation}
\begin{aligned}
v &= \mathbb E^{\mathbb Q}[h(S_{T})] \\
&= \mathbb Q(W_1 < - \frac{b + \mu}{\sigma} = -2) \\
&= \Phi(-2), \\
\end{aligned}
\end{equation}

where $W_1 \sim \mathcal N(0,1)$ and $\Phi(\cdot)$ is the $\mathrm{CDF}$ of the standard Normal Distribution.

Thus, the exace price is shown as follows:



```
import scipy.stats as ss
exact_price = ss.norm(0,1).cdf(-2)
print(exact_price)
```

    0.022750131948179195


- Use $\mathrm{OMC}$ find the price

__Soln:__

By using $\mathrm{OMC}$, we replace $\Phi(-2)$ with $\frac{1}{n} \sum_{i=1}^n I(X_i<-2)$


```
import numpy as np

def OMC(N):         # N is the number of samples
  s = 0
  for i in range(N):
    X = np.random.normal(0,1)
    if (X<-2):
      s += 1
  return s/N
OMC(1000000)
```




    0.022768



- Use $\mathrm{IS}(\alpha)$ find the price
 
 __Soln:__

 Assume $\phi_\alpha(\cdot) \sim \mathcal N(-\alpha,1)$, where $\phi_\alpha(\cdot)$ is the $\mathrm{pdf}$.

Then we can rewrite v which has been discussed in the video:

\begin{equation}
\begin{aligned}
v &= \mathbb E[I(Y<-2) \cdot e^{\frac{1}{2}\alpha^2 + \alpha Y}|Y \sim \phi_\alpha] \\
&\approx \frac{1}{n} \sum_{i=1}^n I(Y_i<-2) \cdot e^{\frac{1}{2}\alpha^2 + \alpha Y_i}, \\
\end{aligned}
\end{equation}

where $Y_i \sim \mathcal N(-\alpha,1).$


```
import numpy as np
import math

def ImpSamp(N,alpha):         # N is the number of samples
  s = 0
  for i in range(N):
    Y = np.random.normal(-alpha,1)
    if (Y<-2):
      s += math.e ** (alpha**2/2 + alpha*Y)
  return s/N
ImpSamp(10000,2)
```




    0.02249178731207442



- Can you show your approach is optimal?



```
err_list = []
for alpha in range(10):
  M = 1000                 # the number of trials to find the minimum of error
  error = []
  for i in range (M):
    error.append(abs(ImpSamp(10000,alpha) - exact_price))     # document every trial's error
  err_list.append(min(error))       # document the minimum of error wrt different alpha
min_err = min(err_list)
print(">>> The minimum error among different α is " + str(min_err) + ", where α = " + str(err_list.index(min_err))) 
```

    >>> The minimum error among different α is 2.0011788371201988e-07, where α = 2


- Prove or demonstrate $\mathrm{IS}$ is more efficient to $\mathrm{OMC}$


```
N = 1000                         # N is the number of samples for both OMC and Importance Sampling
sq_err_omc = 0
sq_err_ImpSamp = 0
K = 10000                         # K is the number of trials for OMC
for i in range(K):
  sq_err_omc += (exact_price - OMC(N))**2
  sq_err_ImpSamp += (exact_price - ImpSamp(N,2))**2
mse_omc = sq_err_omc / K
mse_ImpSamp = sq_err_ImpSamp / K
print(">>> The Mean Square Error of Ordinary Monte Carlo is " + str(mse_omc))
print(">>> The Mean Square Error of Importance Sampling is " + str(mse_ImpSamp))  
```

    >>> The Mean Square Error of Ordinary Monte Carlo is 2.228314286581901e-05
    >>> The Mean Square Error of Importance Sampling is 1.1975332166123286e-06


- After trying two different values of $K$ ($K=1000\ or\ 10000$), the value of exponential will not change. Thus we can say that $\mathrm{IS}$ is more efficient to $\mathrm{OMC}$.

__Numeric Method:__

As we have shown before,
$$ \hat v_{OMC} = \frac{1}{n} \sum_{i=1}^n I(X_i<-2),$$
and
$$ \hat v_{IS} = \frac{1}{n} \sum_{i=1}^n I(Y_i<-2) \cdot e^{\frac{1}{2}\alpha^2 + \alpha Y_i}.$$
Both $v_{OMC}$ and $v_{IS}$ are unbiased, thus we can deduce that

\begin{equation}
\begin{aligned}
\mathrm{MSE}(\hat v_{OMC}) &= Var(\hat v_{OMC}) \\
&= \frac{1}{n} \left( \mathbb E[I(X_i<-2)^2] - (\mathbb E[I(X_i<-2)])^2 \right) \\
&= \frac{1}{n}(\Phi(-2) - \Phi(-2)^2) \\
\mathrm{MSE}(\hat v_{IS}) &=  Var(\hat v_{IS}) \\
&= \frac{1}{n} \left( \mathbb{E} [I(Y_{i} < - 2) e^{\alpha^{2} + 2 \alpha Y_{i}}] - \Phi^{2}(-2) \right) \\
\end{aligned}
\end{equation}

Now we calculate the right hand side of the last equation:

\begin{equation}
\begin{aligned}
\mathbb{E} [I(Y_{i} < - 2) e^{\alpha^{2} + 2 \alpha Y_{i}}] &= \int_{- \infty}^{-2} e^{\alpha^{2}+ 2 \alpha y}  \frac{1}{\sqrt{2 \pi}} e^{- \frac{(y + \alpha)^2}{2}} \mathrm{d} y \\
&= \int_{- \infty}^{-2} \frac{1}{\sqrt{2 \pi}} e^{- \frac{y^{2} - \alpha y - \alpha^{2}}{2}} \mathrm{d} y \\
&= \int_{- \infty}^{-2} \frac{1}{\sqrt{2 \pi}} e^{- \frac{(y - \alpha)^{2}}{2}}  e^{\alpha^{2}} \mathrm{d} y  \\
&= e^{\alpha^{2}} \cdot \Phi(-2-\alpha), \\
\end{aligned}
\end{equation}

Thus, we will have
$$\mathrm{MSE}(\hat v_{IS}) = \frac{1}{n} \left( e^{\alpha^{2}} \Phi(-2-\alpha) - \Phi^{2}(-2) \right).$$

Therefore,
$$\mathrm{MSE}(\hat v_{OMC}) - \mathrm{MSE}(\hat v_{IS}) = \frac{1}{n} \left( \Phi(-2) - e^{\alpha^{2}} \Phi(-2 - \alpha) \right) \ge 0.$$


```
import numpy as np
import scipy.stats as ss

for alpha in range(5):
  diff = ss.norm(0,1).cdf(-2) - np.exp(alpha**2) * ss.norm.cdf(-2-alpha)
  print(diff)
```

    0.0
    0.019080728658526478
    0.021020940734838522
    0.020427370203270685
    0.01398320509620665


8c0be3d53ee1e59832f3adbcf495e031378fccbe
ipynb
Jupyter Notebook
src/Hw11_IS.ipynb
Jun-629/20MA573
addad663d2dede0422ae690e49b230815aea4c70
MIT
src/Hw11_IS.ipynb
Jun-629/20MA573
addad663d2dede0422ae690e49b230815aea4c70
MIT
src/Hw11_IS.ipynb
Jun-629/20MA573
addad663d2dede0422ae690e49b230815aea4c70
MIT
2020-02-05T21:42:08.000Z
2020-02-05T21:42:08.000Z
Qwen/Qwen-72B
1. YES
2. YES