Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

How can I calculate the Total Variation Distance of a transition Matrix? is there any built in function? I've searched all documentation and haven;t found anything.

** More information:

Let me try to explain it better. let's say we have a transition matrix (P), 4x4 that describes the probability of going from a, b, c and d to a, b, c or d in 1 step.

We can calculate the stationary distribution of P, and that's called Pi in the following equation, and Pyx is the probability of going from state y to state x (a,b,c,d):

enter image description here

What I want to do is calculate the Total Variation Distance of P from Pi after n-steps and starting on a given state.

* This is what I have so far:

M = {{0.3, 0, 0.5, 0.2}, {0, 0.4, 0.3, 0.3}, {0.3, 0.2, 0, 0.5}, {0.4, 0.1, 0, 0.5}}

B = Transpose[M]

N[B] // MatrixForm

{eVals, eVecs} = Eigensystem[B]

eVals // MatrixForm

eVecs // MatrixForm

eigenvector = eVecs[[1]]

Print["Stationary Distribution"];

eigenvector/Total[eigenvector]

Print["M after 1 step"];

M2 = MatrixPower[M, 2]

share|improve this question
    
And what would that be in the context of markov chains? I would wild-guess the maximum difference (Max[mat1-mat2]) between the probabilities of events in some time given some initial condition (MatrixPower[Q, n].initialVector) versus in the limit should there be any? –  Rojo Oct 5 '12 at 3:05
1  
Have you tried writing your own function? –  rm -rf Oct 5 '12 at 4:26
    
In the Wikipedia definition, there are two probability distributions P and Q, and the total variation is defined as a function of the two. In your question, what are P and Q? –  bill s Oct 5 '12 at 8:49
    
Let's say that I have a transition matrix (n x n), called P, and that I calculate its stationary distribution, which is an eigenvector with an eigenvalue of 1 (1 x n), and let'sc all it Q. Those are the two inputs of the equation in Wikipedia. –  whynot Oct 5 '12 at 14:11
    
Help me see what I'm missing. The eigenvector has the probabilities of being in a certain state. The matrix the probabilities of "going from one to the other". Those aren't the same events, you can't compare them. That's why the matrix has n^2 elements and the vector n –  Rojo Oct 5 '12 at 14:29

2 Answers 2

up vote 1 down vote accepted

The terms in the total variation distance need to be probability distributions. Here's one way to do it: take any starting vector $m_0$ (with nonzero entries that sum to $1$). You can interpret $m_0$ as the initial distribution of the states. Then $m_1=m_0\cdot P$ is the distribution of the states after one time step, $m_2=m_0\cdot P\cdot P$ is the distribution of the states after two time steps, etc., where $P$ is the transition matrix of the Markov chain. Now you have several distributions and you can apply the total variation (as in the Wikipedia article). This would be

 Total[Abs[m0 - m1]]/2

for the distance from $m_0$ to $m_1$. This is really just the 1-norm, so a more concise expression would be

Norm[m0 - m1, 1]/2

You could also take the distance between any of these distributions and the stationary distribution, as was calculated in stationary distribution of a transition matrix

share|improve this answer
    
Of course, Norm[m0 - m1, 1]/2 is equivalent to ManhattanDistance[m0, m1]/2. –  J. M. Oct 6 '12 at 3:16

You can get the stationary distribution by getting the first eigenvector. Or, if I?m mistaken with this, just use eigensystem and filte by the eigenvector 1.

pi[p_] := Normalize[First@Eigenvectors[p, 1], Norm[#, 1] &]


totVar[m1_, m2_] := Total@Abs[m1- m2]/2

So if p is your transition matrix, m the vector with initial states and n the time

totVarMarkov[p_, m_, n_] := totVar{[pi[p]}, {Flatten@m}.MatrixPower[p, n]]

Code untested

share|improve this answer
    
More compactly: totVar[m1_, m2_] := ManhattanDistance[m1, m2]/2. –  J. M. Oct 6 '12 at 3:15

Your Answer

 
discard

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.