I will publish some Matlab features. I would also very thankful for your comments and recommendations.
First of all I’d like to introduce draft for Markov Chain Builder from sequence of values:
For example we have the following sequence of values (number of jobs per seconds):
{2,4,5,6,7,3,3,3,2,3,4,5,3}
With this scripts you can build transition matrix:
Version 1:
X=unidrnd(4,1,100); % Let’s define the matrix with states
%X=[0 0 1 1 0 2 1 0 1 2 0 1 2 0 1 2 0 1 0 1 5 10]; %or just define sequence manuallyX=X+1;% increment the state because if we have zero state this will make our algorithm wrong
states=unique(X); % states 0 – 1, 1 – 2, 2 – 3
Num_of_states=length(states); %Number of states
F=zeros(Num_of_states,Num_of_states); %Frequency matrix
P=zeros(Num_of_states,Num_of_states); %Prob matrix
row=1;column=1;for i=1:length(X)-1
pos_curr=find(states==X(i)); %use absolute values for states as to algorithm
pos_nxt=find(states==X(i+1));
F(pos_curr,pos_nxt)=F(pos_curr,pos_nxt)+1;
end;
row_sum=sum(F,2); %sum of row in order P<=1
for i=1:length(F)
for j=1:length(F)
if (row_sum(i)==0) continue; % we need this check if sum of rows will be zero so no transition occur and “div by zero”
end;
P(i,j)=F(i,j)/row_sum(i);
end;
end;
%clear i;clear j; clear row_sum; clear pos_curr; clear pos_next;
Version 2: More elegant and was found at (http://www.eng-tips.com/viewthread.cfm?qid=236532&page=3)
BrCo:
This exact question came up on CSSM last year. There are many ways to solve these sorts of problems. My favourite is to use the sparse() function:% Some data with repeated sequences
x=[1 6 1 6 4 4 4 3 1 2 2 3 4 5 4 5 2 6 2 6 2 6];
sparse(x(1:end-1),x(2:end),1)
ans =
(3,1) 1
(6,1) 1
(1,2) 1
(2,2) 1
(5,2) 1
(6,2) 2
(2,3) 1
(4,3) 1
(3,4) 1
(4,4) 2
(5,4) 1
(6,4) 1
(4,5) 2
(1,6) 2
(2,6) 3
Or if you want it in array form, where A(i,j) holds the number of changes from i to j:
A=full(sparse(x(1:end-1),x(2:end),1))
A =
0 1 0 0 0 2
0 1 1 0 0 3
1 0 0 1 0 0
0 0 1 2 2 0
0 1 0 1 0 0
1 2 0 1 0 0

More elegant way was found from mk_stochastic function of BNT toolbox (Written by Kevin Murphy)
A=full(sparse(sequence(1:end-1),sequence(2:end),1));
Z = sum(A,2);
S = Z + (Z==0);
norm = repmat(S, 1, size(A,2));
A = A ./ norm;
A new approach was found in markovFit() function from PMTK toolkit, please refer http://pmtk3.googlecode.com/svn/trunk/docs/synopsis/Markov_models.html
Yet another approach:
A=full(sparse(sequence(1:end-1),sequence(2:end),1));
T=normalize(A,2);