3
$\begingroup$

I need some help to generate a Matlab code in order to do the following question. Can somebody help me in this regard. Any sort of hint that could be helpful will surely be appreciated..

Q: "Simulate $N=25$ trajectories of the Ito Process X satisfying the following SDE

$dXt = \mu X_tdt + \sigma X_tdB_t.$

with $X_0=1$, $\mu=1.5$, $\sigma=1.0$ and their Euler approximations with equidistant time steps of size $\Delta=2^{-4}$ corresponding to the same sample paths of the Wiener process on the time interval $[0,T]$ for $T=1$. Evaluate the absolute error by the statistic defined below

m=$\frac1N$$\sum_{k=0}^N $|$X_{T,k}$-$Y_{T,k}$|

where $X_{T,k}$ and $Y_{T,k}$ respectively are the $k$-th simulated trajectories of Ito process and their Euler approximation corresponding to same sample paths of the Wiener process"

I have created the following code on Matlab for the above question. Can somebody correct me if I'm wrong somewhere.

randn('state',100)

mu=1.5; sigma=1; Xzero=1;

T=1; N=25; dt=T/N;

dW=sqrt(dt)*randn(1,N);

W=cumsum(dW);

Xtrue=Xzero*exp((mu-0.5*sigma^2)*([dt:dt:T])+sigma*W);

Xem=zeros(1,N);

Xem(1)=Xzero+dt*mu*Xzero+mu*Xzero*dW(1);

for j=2:N

Xem(j)=Xem(j-1)+dt*mu*Xem(j-1)+sigma*Xem(j)*dW(j);

end

$\endgroup$
  • $\begingroup$ Your code doesn't make much sense. Is it even runnable? I've submitted an edit to try to clean it up but you should check it. But are $u$ and $mu$ the same? What is the variable lambda? And I'm not sure why you're mixing TeX with Matlab code. Finally, unless you're using a really old version of Matlab, you shouldn't be setting the random number generator's seed with 'state' option. Use the rng function, e.g., rng(100). $\endgroup$ – horchler Jun 29 '14 at 17:30
  • $\begingroup$ Adding to comments by @horchler -- $N=25$ is the number of sample paths not time steps. The time step size is supposed to be $dt = T/K = 1/16 = 2^{-4}.$ Your notation has changed several times. For the Euler approximation you just need to fill a $K \times N$ matrix $X$, where for each path index $n$, $X(k+1,n) = X(k,n) + ...$ according to the algorithm. $\endgroup$ – RRL Jun 29 '14 at 18:46
  • $\begingroup$ @RRL, I dont understand what to do. Im so much confused about this question. $\endgroup$ – user158267 Jun 29 '14 at 18:57
  • $\begingroup$ @user158267: Please format your code properly (four spaces or use the button on the toolbar). In addition to what @RRL said, you're still using mu in place of sigma when you integrate numerically. $\endgroup$ – horchler Jun 29 '14 at 21:37
2
$\begingroup$

The increment of Brownian motion $B_{t+ \Delta }- B_t$ is normally distributed with mean $0$ and standard deviation $\sqrt{\Delta}.$

Generate a sample path using the discrete Euler approximation:

$$X_{k+1}=X_{k} + \mu X_k \Delta + \sigma X_k\sqrt{\Delta}\xi\,\,(k=1,2,...),$$

where $\xi$ is a random number with a standard normal distribution.

To generate random samples for $\xi$, first generate a uniformly distributed random number $r \sim$ U(0,1) and take $\xi = N^{-1}(r)$ where $N$ is the standard normal distribution function.

$\endgroup$
  • $\begingroup$ Could you give me a code to generate this on Matlab? $\endgroup$ – user158267 Jun 29 '14 at 1:42
  • $\begingroup$ Unfortunately I don't have Matlab. I've shown you the algorithm. You can code it in one line. To generate the random numbers you can either follow my suggestion -- or most likely Matlab has a library routine for drawing a normal random number directly. Good luck. Let me know if you have further questions. $\endgroup$ – RRL Jun 29 '14 at 1:57
  • $\begingroup$ I have just added a matlab code in my post. Kindly have a look at it and let me know whether im quite right here or not. $\endgroup$ – user158267 Jun 29 '14 at 14:48
  • $\begingroup$ @horchler, I have just done some editing. Moreover I have read this paper but I dont understand what trajectories mean here. $\endgroup$ – user158267 Jun 29 '14 at 18:15
1
$\begingroup$

Your stochastic differential equation is for geometric Brownian motion. Using your notation, the analytic solution for geometric Brownian motion under Itō's interpretation is

$$X_t = X_0\space\text{exp}\left(\left(\mu-\frac{\sigma^2}{2}\right)t + \sigma B_t\right)$$

This is not what you're using in your Matlab code for Xtrue.

To simulate your system, you can use the Euler-Maruyama for Itō SDEs. @RLL describes this very well so I won't repeat it here. Your Matlab code, however, doesn't appear to implement Euler-Maruyama correctly. There's a mystery parameter lambda, the parameter sigma is missing, and the mu parameter is used in diffusion part rather than the drift. As this appears to be an assignment and I think you're capable fixing this, I'll leave the rest to you.

For further details on SDEs, Brownian motion, and simulating them with Matlab I recommend this excellent paper:

Desmond J. Higham, 2001, An Algorithmic Introduction to Numerical Simulation of Stochastic Differential Equations, SIAM Rev. (Educ. Sect.), 43 525–46. http://dx.doi.org/10.1137/S0036144500378302 http://www.caam.rice.edu/~cox/stoch/dhigham.pdf

The URL to the Matlab files in the paper won't work, use this one: http://personal.strath.ac.uk/d.j.higham/algfiles.html

(Note that the Matlab code used is 13 years old now not meant for efficiency. Some things have changed since, e.g., how to set the random seed.)

$\endgroup$
  • $\begingroup$ Dear horchler can u let me know any of the codes given in the paper by Desmond J.Higham which match match with my case so that I could work on that properly. $\endgroup$ – user158267 Jun 29 '14 at 18:18
  • $\begingroup$ @user158267: The Euler-Maruyama method is discussed in Section 4 of the Higham paper, but perhaps you should read the whole thing and try to understand it rather than taking a shortcut. The paper uses different names for variables so don't blindly copy the code. You need to figure out which correspond to the ones in your SDE. If Matlab programming is partly the issue, perhaps you should step back and spend some time doing a few tutorials. $\endgroup$ – horchler Jun 29 '14 at 21:50
  • $\begingroup$ $N=25$ trajectories just means 25 separate simulation runs (of both the analytical solution, $X_{t,k}$, and the numerical solution calculated via Euler-Maruyama, $Y_{t,k}$), where $k$ is the index of the simulation. You'll want to store the last value of each of these simulations (i.e., $X_{T,k}$ and $Y_{T,k}$) so that you can calculate m. $\endgroup$ – horchler Jun 29 '14 at 21:54
  • $\begingroup$ I have just read the paper of Higham. To my understanding of that paper, I think the code in Listing 6 should be used in my case. Am I right? $\endgroup$ – user158267 Jun 29 '14 at 23:18
  • $\begingroup$ @user158267: Yes, that or Listing 5 looks like where you got your code in the first place. Pretty much none of the parameter/variable names are the same as in your question so you'll need to understand which is which and what each line of the code does to answer your question. If you have programming-related questions –as opposed to math– you can ask at StackOverflow.com/Matlab. But make sure to ask a specific coding question showing what you've tried, what's not working, and what you expect (don't expect people to write the code for you). $\endgroup$ – horchler Jun 30 '14 at 0:22

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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