You can solve DTMC on Octave by using this short code:
P = [
0.6, 0.4, 0, 0;
0, 0.4, 0.6, 0;
0, 0, 0.8, 0.2;
1, 0, 0, 0
]
pis = [P' - eye(size(P)); ones(1, length(P))] \ [zeros(length(P), 1); 1]
Or with SAGE with this code:
P = matrix(RR, 4, [
[0.6, 0.4, 0, 0],
[ 0, 0.4, 0.6, 0],
[ 0, 0, 0.8, 0.2],
[ 1, 0, 0, 0]
])
I = matrix(4, 4, 1); # I; I.parent()
s0, s1, s2, s3 = var('s0, s1, s2, s3')
eqs = vector((s0, s1, s2, s3)) * (P-I); eqs[0]; eqs[1]; eqs[2]; eqs[3]
pis = solve([
eqs[0] == 0,
eqs[1] == 0,
eqs[2] == 0,
eqs[3] == 0,
s0+s1+s2+s3==1], s0, s1, s2, s3)
On both, the result of the steady state probabilities vector is:
pis =
0.245902
0.163934
0.491803
0.098361
I hope it helps.
WBR, Albert.