Looping for Markov Chain

1 visualizzazione (ultimi 30 giorni)
Wendell
Wendell il 25 Ott 2012
I need to run the following:
P = [0.942 0 0 0 0.058 0 0 0; 0 0.982 0 0 0.018 0 0 0; 0 0 0.999 0 0 0.001 0 0; 0 0 0 0.993 0 0.007 0 0; 0 0 0 0 0.987 0 0.013 0; 0 0 0 0 0 0.999 0.001 0; 0 0 0 0 0 0 0.972 0.028; 0 0 0 0 0 0 0 1]; x_0 = [1000; 800; 300; 400; 700; 200; 300; 0]; n = 100; for i = 1:n x = P*x_0; x(i) = P^(i)*x_0; end disp (x(i)) plot (x(i),n)

Risposte (1)

Kye Taylor
Kye Taylor il 25 Ott 2012
Modificato: Kye Taylor il 25 Ott 2012
You have all the pieces... its just awkward to organize (below I use cells) and display.. Try something like
n = 100;
P = cell(1,n);
P{1} = [0.942 0 0 0 0.058 0 0 0; 0 0.982 0 0 0.018 0 0 0; 0 0 0.999 0 0 0.001 0 0; 0 0 0 0.993 0 0.007 0 0; 0 0 0 0 0.987 0 0.013 0; 0 0 0 0 0 0.999 0.001 0; 0 0 0 0 0 0 0.972 0.028; 0 0 0 0 0 0 0 1];
x_0 = [1000; 800; 300; 400; 700; 200; 300; 0];
% I expect x_0 to be a probability
x_0 = x_0/sum(x_0);
% Create P, P^2, P^3,...
for i = 2:n
P{i} = P{1}*P{i-1};
end
g = @(A)A*x_0;
evolvedX = cellfun(g,P,'un',0); % applies P^k to x_0 for k = 1,2,3,...
plot([evolvedX{:}]')
xlabel('Time-step')
ylabel('Probability')
legend('x_0(1)','x_0(2)','x_0(3)','x_0(4)','x_0(5)','x_0(6)','x_0(7)','x_0(8)')
  1 Commento
Wendell
Wendell il 26 Ott 2012
Kye I appreciate very much the help...glad to have someone willing to share their knowledge and skill...

Accedi per commentare.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by