# Can Mathematica plot this financial simulation?

I'd like to run a simple financial simulation with something like the following constraints:

• You start out with $$S$$ dollars. Let's say $$S$$ is $$\1M$$ dollars
• You repeatedly make an investment that costs $$I$$ dollars (let's say $$I = \5K$$) with an $$F$$ percent (let's say $$F = 0.95$$) chance of returning $$\0$$ and an $$(1-F)$$ percent chance of returning $$\%100,000$$ of $$I$$ (or just $$1,000 \cdot I$$).

Is there a way in Mathematica to repeatedly make this investment hundreds of times over time? Could we have Mathematica repeatedly run the simulation and then plot the average of all of its simulation as a line as well?

• Are you sure you didn't mean to write the chances the other way around ie F=.95; also, is there a chance you meant for the returns in the good case to read 1000I instead of 1000S or something similar? – yosimitsu kodanuri May 23 '19 at 17:36
• Yes, you're correct. I adjusted my question above. – George May 23 '19 at 17:39

It can certainly be done using NestList:

s = 100;
i = 0.5;
f = 0.95;
nStep = 100;
sim := NestList[If[RandomReal[] < f, # - i, # - i + 1000 i] &, s,
nStep];


Now sim will generate a simulation of such process at each evaluation. You can use the following to do multiple runs and plot:

runs = Table[sim, {r, 100}];
ListPlot[runs, Joined -> True, PlotRange -> All]


Also plot the means:

ListPlot[Mean @ runs, Joined -> True, PlotRange -> All]

• It's hard to parse the output. Is there a way to adjust the plot to see what's going on? – George May 23 '19 at 17:18
• It looks like the members of sims are growing by 1000x every step. Are you sure this is correct? They should have a 5% chance of growing by 1000x and a 95% of losing the price of the investment (\$5K in this example). – George May 23 '19 at 17:23
• Well that was according to your previous version. Edited. – L.Yu May 23 '19 at 17:52
• Thank you for putting up with my change in question :) – George May 23 '19 at 17:54
• No problem. You can adjust the s,i or f parameters as you wish. – L.Yu May 23 '19 at 17:55

In what follows, the code provided will return a plot of a (user provided) number of random experiments that implement the concept of the investment lottery in the question; the output is a sequence of the amounts outstanding, after a lottery has been paid for and possible profits have been accrued instead of the sequence of winnings.

Define an auxiliary function that generates appropriate random numbers:

randomSequence[distribution_: UniformDistribution[], repetitions_: 10, observations_: 10^6] :=
RandomVariate[distribution, {repetitions, observations}]


Eg. evaluating randomSequence[NormalDistribution[], 2, 10] will return a $$2\times10$$ matrix with random numbers from $$N(0,1)$$.

Define a lottery without using conditionals; make use of the main evaluation loop instead:

(* when rand >= cut-off prob, return a prize *)
lottery[randomValue_, cost_, multiplier_, probability_] /; randomValue >= probability :=
multiplier cost

(* ...otherwise return 0. *)
lottery[__] := 0.


balance represents the amount remaining after committing to participate in an "investment" lottery; so far nothing fancy, but it's possible to modify balance to obtain more "involved" results if needed.

SetAttributes[balance, HoldFirst]
balance[lottery[randomValue_, cost_, multiplier_, probability_], wealth_] :=
wealth - cost + lottery[randomValue, cost, multiplier, probability]


simulation is the main event; it bundles together the simulation results and the display of the results:

SetAttributes[simulation, HoldFirst]
simulation[randomSequence[args___], capital_, investment_, multiplier_, probability_] :=
series = FoldList[
balance[lottery[#2, investment, multiplier, probability], #1] &,
capital, #] & /@ randomSequence[args];

legend = Column[Thread[{"K", "I", "x", "p"} -> {capital, investment, multiplier, probability}]];
opts = {
Frame -> True,
FrameLabel -> {{None, None}, {None, label}}};

Legended[ListLinePlot[series, Apply[Sequence][opts]], legend]
]


Evaluating

 simulation[randomSequence[UniformDistribution[], 10, 10^3], 10^6, 5 10^3, 10^3, .95]


returns