function v = stat(P)
%
% This is a MATLAB function that calculates the stationary probability vector v
% of a Markov chain transition matrix P, i.e., we solve v = v*P .
% We assume the existence of a unique stationary vector.
% For a finite-state Markov chain, the condition is that the chain be irreducible.
% 
% We input the matrix P when we call the function.
% First find the number n of rows in the transition matrix P.
s = size(P);
n = s(1);
%
% There is one redundant equation in the n equations v = vP.
% We fill gap by using the fact that v(1) + ... + v(n) = 1.
% First, we can rewrite v = vP by v(P-I) = z, where I is an identity matrix and z is a vector of zeros.
% We then add a column of 1's to make a new equation
%
I = eye(n); %the identity matrix
z = zeros(1,n); %a row of zeros
w = ones(n,1); %a column of 1's
A = [P-I w];
%
% The desired system of equations is vA = [z 1], where A is n by (n+1)
% We solve it by writing v = [z 1]/A
%
v = [z 1]/A;
%
% Using transposes, we could also write v' = A'\[z 1]'
% We could also use the matrix inverse applied to square matrices.
% That approach is in the other program stationary.m