Mathematica Stack Exchange is a question and answer site for users of Mathematica. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I have a discrete- time and space markov process that I want to evolve given an initial distribution.

To that, I defined

TimeP[initial_, state_, 0] := initial[state];
TimeP[initial_, state_, t_] := TimeP[initial, state, t] = TimeP[initial, state, t - 1] + ...

where the initial distribution would be something like this:

initial[state_] := 0
initial[firstState] := 1

and ... are probabilities to jump from a state to another state.

Now, with the above definition of TimeP, I can correctly evolve the probabilities. However, when I define another initial distribution, say

initial[state_] := 0
initial[firstState1] := 1

The function TimeP gets confused because it only stores the name initial, and not its definition (i.e. even if initial has now a different definition, TimeP returns the computed value from initial's definition of firstState, and not firstState1.

So, my question is: is a way of telling TimeP that its argument's name doesn't matter, only the actual definition that matters? Or, do you know a better way of implementing this?

This is what I would expect: the function initial is defining a distribution, and the distribution's definition is on its probabilities. When I change initial's definition to another values, TimeP would give another result.

(I'm not using Mathematica's DiscreteMarkovProcess because the transition matrix is almost diagonal with most entries equal to 0, and the number of states is large).

share|improve this question
It would seem strange to me to rename your initial state whenever you want to specify a new initial distribution (firstState vs firstState1). Can you not simply assign a different value to initial[firstState]? Also, do you have a value set for firstState? – Jacob Akkerboom Oct 21 '13 at 10:48

It is not a direct answer to your question, but it can helps you.

I think it is more efficient to work with sparse arrays and store the resulting distribution in an array, not in DownValues.

n = 100;

m = SparseArray[{Band[{1, 1}] -> 0.5, Band[{1, 2}] -> 0.25, 
   Band[{2, 1}] -> 0.25}, {n, n}];

init = 50;

p0 = MapAt[1. &, #, init] &@ConstantArray[0., n];

p = NestList[m.# &, p0, 100];

ListPlot3D[p, PlotRange -> {0, 0.2}]

enter image description here

Approach with the memoization:

distr[init_, t_] := distr[init, t] = m.distr[init, t - 1];
distr[init_, 0] := m.MapAt[1. &, #, init] &@ConstantArray[0., n];

ListLinePlot[{distr[50, 100], distr[75, 100], distr[50, 200]}]

enter image description here

share|improve this answer

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

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