function [chain,state] = markovsimul(PI,s,n,s0,seed)
% Simula uma cadeia de Markov
%  [chain,state]=markovsimul(PI,s,n,s0,seed)
% INPUTS:
%   PI	: Matriz de Transição
%   s	: Vetor de Estados
%   n	: Tamanho da simulação
%   s0	: indice do estado inicial
%   seed: semente para geração de numeros aleatórios
%
% OUTPUTS:
%   chain : values da cadeia de markov simulada
%   state : indices dos estados correspondentes a cadeia de markov simulada

[rpi,cpi]=size(PI);
s=s(:);
cum_PI = [zeros(rpi,1) cumsum(PI,2)]; % probabilidade acumulada em cada uma das linhas

rand('state',seed) % define seed do gerador de numeros aleatorios para permitir replicabilidade
sim		= rand(n,1); % sorteia n numeros aleatorios de uma uniforme 0-1
state		= zeros(n,1); % vetor que guardara os estados da cadeia simulada
state(1) 	= s0; % estado inicial dado


for k=2:n
    state(k) = find(sim(k)>=cum_PI(state(k-1),:),1,'last');
end
chain=s(state);