I want to price exotic options under the exponential VG model and Merton's model to compare both models.

To price exotics under Merton's model, I have written the code below. The output is the price of a Call option, Asian, Digital and Up and in Barrier Call option. However, the use of loops leads to a very slow computation. Is there a clever way to not use loops here? In case of the VG model, I can do it but in this case, I do not see it.

function [Call,Asian,Digital,UIBP] = ExoticPricingMerton(S0,K,mu,delta,lambda,sigma,r,q,Maturity,H)

ht = 1/252; %trading days
P = 10^3; %Number of simulations
grid = (0:ht:Maturity); 
N = length(grid);

omega = r-q-((1/2)*sigma^2+lambda*(exp(mu+(1/2)*delta^2)-1)); 

S = zeros(P,N);
S(:,1) = S0; 

for i=1:P
    for j=2:length(grid)
        N = poissrnd(lambda*ht);
        J = cumsum([0, normrnd(mu,delta,1,N)]); 
        Z = normrnd(0,1);
        S(i,j) = S(i,j-1)*exp(omega*ht + sigma*sqrt(ht)*Z + J(end));

%European Call option
A = max(S(:,end)-K,0); 
Call = exp(-r*Maturity)*(1/P)*sum(A); 

%Asian option
A = max(mean(S,2) - K,0);
Asian = exp(-r*Maturity)*(1/P)*sum(A);

%Digital price
A = max(S(:,end) - K, 0)./(S(:,end)-K);
Digital = exp(-r*Maturity)*(1/P)*sum(A);

%Up-and-in out Barrier
A1 = (max(S,[],2)-H)./abs(max(S,[],2)-H);
A2 = max(A1,0);
A = (max(S(:,end)-K,0)).*A2;
UIBP = exp(-r*Maturity)*(1/P)*sum(A); 


  • $\begingroup$ Hi, I have been solving this issue recently. What caused troubles in my program was the random draws as they turned out to be very slow. I solved it generating random matrix rand(m,n) where m is number of paths and n is number of grid points. For the poisson process u can use poissrnd(lambda,m,n). Then in forloops just select whatever discretization scheme you want. $\endgroup$ Commented Aug 14, 2016 at 17:31
  • $\begingroup$ Hi, thanks for your answer. I am aware that one can generate random numbers in matrix form. However, I'm not sure how to do it for compound Poisson random variables. $\endgroup$
    – user39039
    Commented Aug 14, 2016 at 17:42
  • $\begingroup$ the "J = cumsum([0, normrnd(mu,delta,1,N)]);" does it sum over N draws? You don't seem to save the draws so why not J = sqrt(N) * normrnd(mu,delta,1,1)? $\endgroup$
    – Mats Lind
    Commented Aug 14, 2016 at 17:55
  • $\begingroup$ @MLind the cumsum does indeed sum over the N draws. Thanks for noticing, that is indeed simpler. However, I do not think it will reduce the overall computation time. $\endgroup$
    – user39039
    Commented Aug 14, 2016 at 18:03
  • 1
    $\begingroup$ @user39039 Generally, Monte-Carlo Method is time-consuming. $\endgroup$
    – user16651
    Commented Aug 14, 2016 at 19:59

2 Answers 2


You're using a wrong tool for the job. Write your Monte Carlo in a faster language (Java would probably suffice, if not than C++ which is standard for such things). Then you will be able to efficiently generate more than 1000 paths. In fact, doing Monte Carlo derivatives pricing with 1000 paths is worthless. Your results are, most probably, very inaccurate. Read a good book on Monte Carlo pricing before venturing further and wasting your time.

  • $\begingroup$ I agree with your sentiment, but are there any connectors between Matlab and C++/Java? $\endgroup$
    – Owe Jessen
    Commented Dec 30, 2016 at 15:16
  • $\begingroup$ @owejessen yes, you'll likely need to create the wrapper though, but you can for sure do it. $\endgroup$
    – will
    Commented Dec 30, 2016 at 16:20
  • $\begingroup$ @will That's good to know, not for myself, but for OP. I live in the R world, and tried to replicate the functions in R / RCPP. It certainly takes a lot of time to run in pure R (as in Matlab), but I'm too dumb to recode in cpp. $\endgroup$
    – Owe Jessen
    Commented Dec 30, 2016 at 16:56
  • $\begingroup$ If you want to use cpp in R, you can jsut write a COM interface to the cpp code and call it from R. Calling cpp from MatLab, from R. You can also create .Net wrappers using CLI - we use this to call from C# inside an excel-DNA addin. We've used all of these approaches, and they all work nicely. $\endgroup$
    – will
    Commented Dec 30, 2016 at 17:02
  • $\begingroup$ It's fairly easy to call Java code from Matlab. $\endgroup$
    – quant_dev
    Commented Dec 31, 2016 at 16:25

I assume you are using MatLab.

You may consider pre-generating all 1,000 random numbers once before for-loop by exploiting array coding.

Another approach, have you ever tried using Quasi Monte Carlo?

Generating Quasi-Random Numbers

QMC ensures faster convergence and MatLab has functions that can generate quasi-random sequence very fast (a billion under a second).

  • $\begingroup$ For 1000 paths, I don't think the difference between pseudo- and quasi-random numbers will be noticeable. But I might be wrong. $\endgroup$
    – quant_dev
    Commented Feb 2, 2017 at 19:48
  • $\begingroup$ @quant_dev You are right. Quasi will be slower, but it will guarantee better convergence. And I don't think 1000-path will be enough for practical use anyways. $\endgroup$
    – Woraphon T
    Commented Feb 3, 2017 at 0:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.