17
$\begingroup$

What is the canonical way of simulating discrete time stochastic dynamical systems in Mathematica using the new functionality of Random processes?

To take a concrete example, lets consider the optimal gambling problem. A gambler comes to a casino with an initial fortune $x_1$ and let $X_n$ denote his fortune at time $n$. At each time he bets a fraction $\alpha x_n$. At each time he wins with probability $p$ and loses with probability $(1-p)$. Thus, the dynamics of the system can be written as:

$$ X_{n+1} = \begin{cases} (1+\alpha) X_n, & \text{with probability } p \\ (1-\alpha) X_n, & \text{with probability } 1-p \\ \end{cases} $$

I want to simulate this system for 20 time steps with $x_1 = 12$, $p=0.8$, and then plot 50 sample paths and the mean value of $X_{20}$.

|improve this question
$\endgroup$
15
$\begingroup$

I would proceed like the following. It will be natural to propose that the win-event occurs following BinomialDistribution with probability $p=0.8$ so that we can use the built-in BinomialProcess to simulate the win and losses in $20$ time steps for $50$ sample paths.

timstep = 20;
win = BinomialProcess[.8];
samplepaths=50;
process = RandomFunction[Evaluate[win], {0, timstep - 1}, samplepaths];
ListLinePlot[process, InterpolationOrder -> 0, Frame -> True,
PlotStyle -> {{Dashed, Opacity[.9]}}]

enter image description here

Then just use FoldList on the TemporalData to get your result simulation of $50$ sample paths and the mean in a simple functional style!

x1 = 12;
alpha = .15; (* assume some value for alpha *)
WinEvents = process["States"];
ST[WinEvents_?ListQ, x1_, alpha_, SimTime_] := Module[{},
   FoldList[Max[If[#2 == 1, (1 + alpha) #1, (1 - alpha) #1], 0.] &,
            x1,
            Prepend[Differences@WinEvents, First@WinEvents]]
  ];
simdat = ST[#, x1, alpha, timstep] & /@ WinEvents;

For visualization

Show[ListLinePlot[Legended[Mean[simdat], "Mean"], PlotRange -> All, 
PlotStyle -> {{AbsoluteThickness[3.5], Opacity[.9], Red}}], 
ListLinePlot[simdat, PlotRange -> All, 
PlotStyle -> {{Dashed, Opacity[.7]}}], 
ListPlot[simdat, PlotRange -> All], Frame -> True]

enter image description here

|improve this answer
$\endgroup$
  • $\begingroup$ Thanks. this is exactly what I was looking for. I did not know about process["States"] and now I can read up more about it. $\endgroup$ – Aditya Apr 15 '13 at 9:57
  • $\begingroup$ @Aditya You are welcome! Mathematica doc has got lots of hidden gems..wish you a happy treasure hunting ;) $\endgroup$ – PlatoManiac Apr 15 '13 at 10:01
  • $\begingroup$ An alternative option is to use BernoulliProcess instead of BinomialProcess. Then, one would not need to take the differences in ST function. $\endgroup$ – Aditya Apr 15 '13 at 10:11
  • $\begingroup$ You are right! I randomly chose a distribution for this example.. $\endgroup$ – PlatoManiac Apr 15 '13 at 11:10

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.