Monte Carlo Pricing of Options Across Multiple Underlyings

# Montecarlo pricing


# Montecarlo pricing












I have some problems with the Montecarlo simulation to price a generic Call option. I want to explain something regarding MC simulation with a simple cases, and after that I am going to talk about my problem.

- Montecarlo - Simple case: considering a set of parameters: S=1, K=1, sigma=0.5, r=0, T=1, N=1000 (number simulation MC); the montecarlo with this example works in this way:

1.1 $S$ = [1,1, ..., 1] I am going to repeat the value of underlying a number of times equal to $N$

1.2 $X$ = $S e^{Z}$ where $Z$ = ($1\times N$) vector of Browniam motion --> I get a vector ($1\times N$) where I multiply each value of $S$ with each value of $e^{Z}$

1.3 In the vector $X$, I take the $\max(X-K,0)$ so between each value inside the vector minus K, and zero

1.4 Finally, I find the payoff, that is the average of the vector.

- Montecarlo - Another case: considering a set of parameters: S=[1.1,1.2], K=1, sigma=0.5, r=0, T=1, N=1000 (number simulation MC); the montecarlo with this example works in this way:

2.1 $S = \begin{pmatrix} 1.1 & \dots & 1.1\\ 1.2 & \dots & 1.2 \end{pmatrix}$ ($2 \times N$); so we repeat the value of underlying a number of times equal to $N$

1.2 $X$ = $S e^{Z} $ = $\begin{pmatrix} 1.1 \exp{Z_1} & \dots & 1.1\exp{Z_N}\\ 1.2 \exp{Z_1}& \dots & 1.2\exp{Z_N} \end{pmatrix} $ where $Z$ = ($1\times N$) vector of Browniam motion --> I get a vector ($2\times N$) where I multiply each value of first row $S$ with each value of $e^{Z}$, the same for the second row of $S$

1.3 In the vector $X$, I take the $\max(X-K,0)$

1.4 Finally, I find the two payoffs, that is the average of the first row (for first payoff) and average of second payoff (for second row)

Now i can explain my problem: How can I find the payoff, if I have both $S$ and $\sigma$ that are vectors? for example, $S$=[1.1,1.2], $\sigma$=[0.5,0.6]

I have tried in this way, but I think it is wrong..

- $S = \begin{pmatrix} 1.1 & \dots & 1.1\\ 1.2 & \dots & 1.2 \end{pmatrix}$($2 \times N$)

- $\sigma = \begin{pmatrix} 0.5 & \dots & 0.5\\ 0.6 & \dots & 0.6 \end{pmatrix}$($2 \times N$)

- I generate two Brownian motion (because I have two values of sigma) = $Z = \begin{pmatrix} BM_{11} & \dots & BM_{1N}\\ BM_{21} & \dots & BM_{2N} \end{pmatrix}$

- $X = S e^{Z} = \begin{pmatrix} \begin{pmatrix} 1.1 e^{BM_{11}} & \dots & 1.1 e^{BM_{1N}}\\ 1.2 e^{BM_{11}} & \dots & 1.2 e^{BM_{1N}} \end{pmatrix} \\ \begin{pmatrix} 1.1 e^{BM_{21}} & \dots & 1.1 e^{BM_{2N}}\\ 1.2 e^{BM_{21}} & \dots & 1.2 e^{BM_{2N}} \end{pmatrix} \end{pmatrix}$

- After that, I take the maximum as before, and the average, but in this case I obtain 4 payoff! And for this reason I am not sure about this method..

## Answer by Kermittfrog (score 5, accepted)

https://quant.stackexchange.com/a/66576

Commonly, you do not use 'pure' matrix algebra when formulating a Monte Carlo valuation setup.

If your options are of the European type, and you truly want to price all options in one go, you could go as follows. Let $M$ denote the number of simulations $m=1\ldots M$, let $N$ denote the number of underlyings $n=1\ldots N$, and let $K_m$ denote each strike.

Let $\Sigma$ denote the covariance matrix of your asset returns, i.e. $\Sigma_{i,j}=\sigma_i\sigma_j\rho_{i,j}$ and $\Sigma_{ii}=\sigma_i^2$. Then, $C$ denotes the Cholesky decomposition of $\Sigma$, i.e. $CC^T=\Sigma$. Note that $\Sigma$ is $N \times N$ and $C$ is $N \times N$ as well.

To stick with your example, let $\sigma\equiv\begin{pmatrix}\sigma_1\\\sigma_2\\\ldots\\\sigma_N\end{pmatrix}$

If you have no correlation between assets, then $$ \begin{align} \Sigma&=\begin{pmatrix}\sigma_1^2&0&\ldots&0\\ 0&\sigma_2^2&\ldots&0\\0&0&\ldots&0\\0&0&\ldots&\sigma_N^2\\ \end{pmatrix}=\mathrm{Diag}\left(\sigma_1^2,\sigma_2^2,\ldots,\sigma_N^2\right)=\mathrm{Diag}\left(\sigma\right)\mathrm{Diag}\left(\sigma\right) \end{align} $$ and, of course, $$C=\begin{pmatrix}\sigma_1&0&\ldots&0\\ 0&\sigma_2&\ldots&0\\0&0&\ldots&0\\0&0&\ldots&\sigma_N\\ \end{pmatrix}=\mathrm{Diag}(\sigma)$$

Disregarding interest rates and dividend yields, let's introduce the $1\times M$ vector $e$ which consists of ones, only. Given an $N\times M$ matrix of standard normal variates $U$, we can now set $$ X=-\frac{1}{2}\sigma\otimes\sigma e T+ \sqrt{T}CU $$ where $\otimes$ denotes element-wise multiplication, and $T$ denotes the time to expiry of your option. The dimension of $X$ is then $N \times M$ as well.

Ultimately, we find

$$ S_T=S_0\otimes e^{X} $$

and $$ P=\max(S_T-Ke,0) $$ with $S_0$ the $N\times 1$ vector of initial prices and the exponential operator $e$ to be understood in an element-by-element way. Furthermore, $\max$ is to be understood in an element-by-element way as well, and $K$ is the vector of strikes. You can then average over each row of the $N\times M$ option payoff matrix $P$.

Again, your programming language of choice should have some helpers and common interpretation of scalar-by-vector, scalar-by-matrix and vector-by-matrix multiplication; and you will most probably not follow this setup in totality.