2
$\begingroup$

I am trying to construct a recursion, g[t], where

  • g[t] = a[t] if g[t-1] = 1
  • g[t] = b[t] if g[t-1] = 2.

Additionally a and b are defined as:

a[t_] := If[0 <= RandomReal[] <= 0.94, 1, 2]
b[t_] := If[0 <= RandomReal[] <= 0.9, 2, 1]

So then I want (I've also tried it where I replace a and b with their original defintions):

g[t_]:= Piecewise[{{a[t], g[t - 1] = 1}, {b[t], g[t - 1] = 2}}]
g[1] = 1

This, however, is not working.

$\endgroup$
4
  • 3
    $\begingroup$ Welcome here. You are aware that g[t-1]=1 does not mean is equal to? Therefore, you first two bullet points don't make sense to me, because after the if I expect a true/false. $\endgroup$
    – halirutan
    Commented Sep 30, 2013 at 12:27
  • $\begingroup$ My first two bullet points are not actually written into my code; I was just trying to give a general idea as to what I'm trying to do. The parts in boxes are what I'm actually working with. But do you mean that my conditions, g[t-1] = 1 or 2, defined in my piecewise function do not make sense? I've noticed that if i choose a simple condition, like t > 1, then my piecewise function works. $\endgroup$ Commented Sep 30, 2013 at 15:36
  • $\begingroup$ Possible duplicate: 18697. $\endgroup$ Commented Oct 1, 2013 at 15:48
  • $\begingroup$ In my opinion, that one is indeed duplicate. However, maybe it is not too localized, as I guess this is a great place to learn for the first time the difference between == and =. I suggest we close this Q&A (but I'm still glad I got to write an answer!) $\endgroup$ Commented Oct 1, 2013 at 15:50

2 Answers 2

2
$\begingroup$

As halirutan has pointed out, the main problem with your code is the equals sign that in the second elements (conditions) of the lists in the list in the Piecewise function. This should be ==.

The t in a[t] and b[t] suggests to me that you want to store values of these random sequences, rather than generating them over and over.

Also you can use RandomChoice, rather than your construction with RandomReal.

To store your values in the code below, I have used something that looks like a technique called memoization. For a and b, this changes the outcome, as things are generated only once, so it is not really memoization. But given a and b, calculating g from its stored previous value rather than going all the way back to t=1 is proper memoization. Note that in g[t_]:=g[t]=body the brackets should be read to be g[t_]:= (g[t]=body). Search this site for memoization if you are not sure how it works :).

The code then becomes

ClearAll[a, b, g]
a[t_] := a[t] = RandomChoice[{0.94, 0.06} -> {1, 2}]
b[t_] := b[t] = RandomChoice[{0.9, 0.1} -> {1, 2}]
g[t_] := g[t] = 
  Piecewise[{{a[t], g[t - 1] == 1}, {b[t], g[t - 1] == 2}}]
g[1] = 1;

Note that more improvements can be made. For example, it is much faster to first generate say 100000 random variables all at once, than to generate them on the fly. But maybe that is something for a future question.

$\endgroup$
1
  • $\begingroup$ I did not search for duplicates, which I am sure there are. Do not hesitate to close this Q&A or whatever, do not let the presence of my answer or any upvotes for it discourage you to close/delete :P. $\endgroup$ Commented Sep 30, 2013 at 16:44
1
$\begingroup$

I originally posted an incorrect post. This may have been better as a comment and does not address the code correct. However, I note that if the aim is a simulation, this two state discrete Markov process staring in state 1 could be coded using Mathematica builtin functions:

markov = DiscreteMarkovProcess[{1, 0}, {{0.94, 0.06}, {0.1, 0.9}}]

Underlying graph:

el = EdgeList[Graph[markov]];
prop = PropertyValue[{Graph[markov], #}, "Probability"] & /@ el;
Graph[markov, EdgeLabels -> Thread[el -> prop], EdgeLabelStyle -> 12]

enter image description here

Random function (e.g 10 samples):

ranf=RandomFunction[markov,{0,10}]

This creates a TemporalData object

Do probability, probability in state 2 after 100 time steps

 Probability[x[100] == 1, x \[Distributed] markov]

(=0.625)

Compute Stationary Distribution:

stat = StationaryDistribution[markov];
$\endgroup$

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.