Simulate a five-state absorbing Markov chain

I am trying to solve this problem :

this my code so far and getting errors:

N=1;      % Number of experiments
n=15;       % Number of transitions to be computed
X=char(zeros(n,N));     % Each column of X represents one of the N experiments
S=char(n,1);            % Initialize the state array
M=zeros(n,5);           % M contains the experimental probabilities for states S & R & N
p00=1;     p01=0;    p02=0;      p03=0;           p04=0;
p10=0.3;   p11=0;    p12=0.7;    p13=0;           p14=0;
p20=0;     p21=0.5;  p22=0;      p23=0.5;         p24=0;
p30=0;     p31=0;    p32=0.6;    p33=0;           p34=0.4;
p40=0;     p41=0;    p42=0;      p43=0;           p44=1;
for j=1:N
s0=randi(4,1);
S(1)=s0;  
for k=1:n-1
r=rand(); s=S(k);
if s=='1'
  if r<=p12, S(k+1)='2'; end 
elseif s=='2'
  if r<=p21, S(k+1)='1'; elseif r>p23, S(k+1)='3'; end
elseif s=='3'
  if r<=p32, S(k+1)='2'; elseif r>p34, S(k+1)='4'; end
elseif s=='4'
  if r<=p44, S(k+1)='4';end
end
end
X(:,j)=S;
end

    % code
  %end
for j=1:n
    x=X(j,:);
    ma=length(find(x==1));
    mb=length(find(x==2));
    mc=length(find(x==3));
    md=length(find(x==4));
    M(j,:)=[ ma mb mc md ];
end
%
nv=0:n-1;
figure(1);
plot(nv, M,'*:');
title('Simulation results -- States A & B');
xlabel('Step number');
ylabel('State');
%

Answer 1

[Akira Agata]

If you have Statistics and Machine Learning Toolbox, you can do this much easier, like:

% Transition matrix
trans = [...
1   0   0   0   0;...
0.3 0   0.7 0   0;...
0   0.5 0   0.5 0;...
0   0   0.6 0   0.4;...
0   0   0   0   1];

% To set the initial state to '2'
trans_hat = [...
0 0 0 1 0 0;
zeros(5,1) trans];

emis = ones(6)/6;
[~,states] = hmmgenerate(15,trans_hat,emis,...
'Statenames',{'tmp','0','1','2','3','4'});

% Show the result
figure
plot(str2double(states),'o-')
ylim([-0.2 4.2])
yticks(0:4)
xticks(0:15)

[MCT001]

I can't use the Statistics and Machine Learning Toolbox:(

[MCT001]

and what about if I want the state to be random? how can I change that ?

[MCT001]

Please see the question part 2

[Akira Agata]

OK. Then, how about the following?

trans = [...
1   0   0   0   0;...
0.3 0   0.7 0   0;...
0   0.5 0   0.5 0;...
0   0   0.6 0   0.4;...
0   0   0   0   1];

nTrans    = 15;             % Number of transisiton
histState = zeros(1,nTrans);% History of visited state
initState = randi(5);       % Randomly select the initial state
currState = initState;

for kk = 1:nTrans
currState = find(rand() <= cumsum(trans(currState,:)),1);
histState(kk) = currState;
end

figure
plot(0:15,[initState,histState]-1,'o-') % To adjust the state label to 0~4
xlabel('# of transision','FontSize',12)
ylabel('State','FontSize',12)
ylim([-0.2 4.2])
yticks(0:4)
xticks(0:15)