I hope this is the right section for this kind of questions.

I am trying to simulate, with MATLAB, a diffusion model starting from a Random Walk. I am using a Random Walk with information increment X normally distributed ($\mu, \sigma$ ). I also have a boundary $\alpha $, and $\alpha > \mu$. The starting point is 0.

If I understood this right, this should be an approximation of the Wiener Process. As wikipedia says (http://en.wikipedia.org/wiki/Inverse_Gaussian_distribution), the first passage time for a fixed level $\alpha > 0 $by $X_t$ is distributed according to an inverse-gaussian: $ T_\alpha = \inf\{ 0 < t \mid X_t=\alpha \} \sim IG(\tfrac\alpha\nu, \tfrac {\alpha^2} {\sigma^2}).\, $

What I am trying to do is to simulate a Random Walk and to get the first passage time distribution, verifying that it is actually a Inv. Gaussian with those parameters.

This is the code I have done: http://pastebin.com/E1N58sJ4

Notice that the myHist function is commented, but I normally use it: in that function I fit the resulting distribution to an inverse gaussian. Then I compare the fitted parameters with the two parameters obtained by the formula showed above, $\mu= \alpha/v $ and $ \lambda= \alpha^2/\sigma^2 $

However, the two results are NOT THE SAME, and they differ consistently across simulation. For example, with the parametrers used in the code:

Simulation - mu:29.1771 s: 116.7757 Expected - mu:26.3158 s: 100

This difference is consistent across repetition of simulation.

Can anyone spot the mistake?


closed as off-topic by Peter Flom Mar 28 '14 at 11:32

  • This question does not appear to be about statistics within the scope defined in the help center.
If this question can be reworded to fit the rules in the help center, please edit the question.

  • $\begingroup$ This question appears to be off-topic because it was already asked on the math site $\endgroup$ – Peter Flom Mar 28 '14 at 11:32
  • $\begingroup$ I already flagged it to a moderator. I first opened it here, and then I realized it was better suited for math exchange. $\endgroup$ – Vaaal Mar 28 '14 at 16:51

Your code uses a quite low sampling rate ($\Delta t = 1$) for simulating the Wiener process. Consider a case where the process crosses the boundary but then returns back to below $\alpha$ before the next sample. This causes $T_\alpha$ to be overestimated, as the first noticed boundary-crossing will occur later than the actual first boundary-crossing.

Indeed, by modifying your code to use $\Delta t = 0.01$, I got a mean of simulated $T_\alpha$s much closer to the analytically obtained 26.3. (I did not check the other parameter as you did not describe how you fit the distribution).

  • $\begingroup$ Hi Juho, thank you for your answer. How would you actually implement that? By the way, I fit the distribution with a Max. Lik. Est. $\endgroup$ – Vaaal Mar 27 '14 at 18:06

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