MC2STAT Calculates the stationary distribution for a Markov chain. CALL: ro = mc2stat(P); ro = Stationary distribution. [1xr] P = transition matrix. [rxr] Example: F = magic(5) P = mat2tmat(F) ro = mc2stat(P) See also mat2tmat, mc2reverse, smc2stat, mctp2stat
Display message and abort function. | |
Convert number to string. (Fast version) | |
Display warning message; disable or enable warning messages. |
Script to computer exercises 2 | |
Script to computer exercises 4 | |
Simulates a Markov chain. | |
Calculates asymmetric rainflow matrix for a MCTP. | |
Calculates the stationary distribution for a MCTP. | |
Demo for switching AR(1)-processes. | |
Rainflow matrix for Switching Markov Chains of Turning Points. | |
Calculates the rainflow matrix/intensity for a switching Markov chain. | |
Simulates a Switching Markov chain with state space. | |
Calculates the asymmetric rainflow matrix for a SMCTP. | |
Calculates the rainflow matrix for a SMCTP. | |
Stationary distributions for a switching MCTP. | |
Simulates a switching Markov chain of turning points, | |
Quick test of the routines in module 'markov' |
001 function [ro,PP]=mc2stat(P) 002 % MC2STAT Calculates the stationary distribution for a Markov chain. 003 % 004 % CALL: ro = mc2stat(P); 005 % 006 % ro = Stationary distribution. [1xr] 007 % 008 % P = transition matrix. [rxr] 009 % 010 % Example: 011 % F = magic(5) 012 % P = mat2tmat(F) 013 % ro = mc2stat(P) 014 % 015 % See also mat2tmat, mc2reverse, smc2stat, mctp2stat 016 017 % Tested on Matlab 5.3 018 % 019 % History: 020 % Revised by PJ 23-Nov-1999 021 % updated for WAFO 022 % Created by PJ (Pär Johannesson) 1998 023 % from 'Toolbox: Rainflow Cycles for Switching Processes V.1.0' 024 % previous name: 'statf_mc' 025 026 % Check input arguments 027 028 ni = nargin; 029 no = nargout; 030 error(nargchk(1,1,ni)); 031 032 n=length(P); % number of states 033 034 PP=P; 035 036 % Set negative elements to zero 037 [I,J] = find(PP<0); 038 if length(I) ~= 0 039 warning(['Warning: Negative elements in P. Setting to zero!']); 040 for k = 1:length(I) 041 PP(I(k),J(k)) = 0; 042 end 043 end 044 045 J =[]; 046 NoDeleted = 0; NoDeleted0 = -1; 047 while NoDeleted ~= NoDeleted0 048 049 % Check that the rowsums of PP are equal to 1 050 051 sumPP = sum(PP'); 052 if sum(abs(sumPP-1) <= n*eps) ~= length(PP) % Number of rowsums equal to 1 053 warning(['Warning: Rowsums of P not equal to 1. Renormalizing!']); 054 for i = 1:length(PP) 055 if sumPP(i) == 0 056 if isempty( find(J==i) ) % Not allready deleted 057 warning(['Warning: Rowsum ' num2str(i) ' of P equals zero. Deleting state ' num2str(i) '!']); 058 J = [J i]; 059 PP(:,i) = zeros(n,1); 060 end 061 else 062 PP(i,:) = PP(i,:)/sumPP(i); 063 end 064 end 065 end 066 067 NoDeleted0 = NoDeleted; 068 NoDeleted = length(J); 069 070 end 071 072 I=1:n; I(J)=[]; 073 P1 = PP(I,I); % Delete states 074 075 076 % Compute stationary distribution 077 078 r=size(P1,1); 079 P1=P1'; 080 M=[eye(r-1) zeros(r-1,1)]; 081 P1=P1(1:r-1,:); 082 P1=P1-M; 083 P1=[P1 ; ones(1,r)]; 084 a=[zeros(r-1,1); 1]; 085 roP1=P1\a; 086 %roP1=pinv(P1)*a; 087 roP1=roP1'; 088 089 % Set ro(i)=0 for 'singular states' 090 091 ro = zeros(1,n); 092 ro(I)=roP1; 093 094
Comments or corrections to the WAFO group