Mercurial > forge
changeset 9733:1a4d4ff140c7 octave-forge
bug fixes
author | mmarzolla |
---|---|
date | Fri, 16 Mar 2012 14:08:44 +0000 |
parents | 0f2ce9bc2525 |
children | d181b0e4c274 |
files | main/queueing/inst/dtmc.m main/queueing/inst/dtmc_check_P.m main/queueing/inst/dtmc_fpt.m main/queueing/inst/dtmc_is_irreducible.m main/queueing/inst/dtmc_mtta.m |
diffstat | 5 files changed, 60 insertions(+), 32 deletions(-) [+] |
line wrap: on
line diff
--- a/main/queueing/inst/dtmc.m Fri Mar 16 11:24:04 2012 +0000 +++ b/main/queueing/inst/dtmc.m Fri Mar 16 14:08:44 2012 +0000 @@ -24,14 +24,16 @@ ## @cindex Discrete time Markov chain ## @cindex Markov chain, stationary probabilities ## @cindex Stationary probabilities +## @cindex Markov chain, transient probabilities +## @cindex Transient probabilities ## -## With a single argument, compute the steady-state probability vector -## @code{@var{p}(1), @dots{}, @var{p}(N)} for a -## Discrete-Time Markov Chain given the @math{N \times N} transition -## probability matrix @var{P}. With three arguments, compute the -## probability vector @code{@var{p}(1), @dots{}, @var{p}(N)} -## after @var{n} steps, given initial probability vector @var{p0} at -## time 0. +## Compute steady-state or transient state occupancy probabilities for a +## Discrete-Time Markov Chain. With a single argument, compute the +## steady-state occupancy probability vector @code{@var{p}(1), @dots{}, +## @var{p}(N)} given the @math{N \times N} transition probability matrix +## @var{P}. With three arguments, compute the state occupancy +## probabilities @code{@var{p}(1), @dots{}, @var{p}(N)} after @var{n} +## steps, given initial occupancy probability vector @var{p0}. ## ## @strong{INPUTS} ## @@ -145,6 +147,15 @@ %! p = dtmc(P, 100, [1 0]); %! assert( plim, p, 1e-5 ); +%!test +%! P = [0 1 0 0 0; \ +%! .25 0 .75 0 0; \ +%! 0 .5 0 .5 0; \ +%! 0 0 .75 0 .25; \ +%! 0 0 0 1 0 ]; +%! p = dtmc(P); +%! assert( p, [.0625 .25 .375 .25 .0625], 10*eps ); + %!demo %! a = 0.2; %! b = 0.15;
--- a/main/queueing/inst/dtmc_check_P.m Fri Mar 16 11:24:04 2012 +0000 +++ b/main/queueing/inst/dtmc_check_P.m Fri Mar 16 14:08:44 2012 +0000 @@ -17,14 +17,14 @@ ## -*- texinfo -*- ## -## @deftypefn {Function File} {[@var{result} @var{err}] =} dtmc_check_P (@var{P}) +## @deftypefn {Function File} {[@var{r} @var{err}] =} dtmc_check_P (@var{P}) ## ## @cindex Markov chain, discrete time ## -## If @var{P} is a valid transition probability matrix, return -## the size (number of rows or columns) of @var{P}. If @var{P} is not -## a transition probability matrix, set @var{result} to zero, and -## @var{err} to an appropriate error string. +## Check if @var{P} is a valid transition probability matrix. If @var{P} +## is valid, @var{r} is the size (number of rows or columns) of @var{P}. +## If @var{P} is not a transition probability matrix, @var{r} is set to +## zero, and @var{err} to an appropriate error string. ## ## @end deftypefn @@ -40,19 +40,19 @@ endif result = 0; + err = ""; if ( !issquare(P) ) err = "P is not a square matrix"; return; endif - if ( any(any(P <0)) || norm( sum(P,2) - 1, "inf" ) > epsilon ) + if ( any(any(P<-epsilon)) || norm( sum(P,2) - 1, "inf" ) > epsilon ) err = "P is not a stochastic matrix"; return; endif result = rows(P); - err = ""; endfunction %!test %! [r err] = dtmc_check_P( [1 1 1; 1 1 1] );
--- a/main/queueing/inst/dtmc_fpt.m Fri Mar 16 11:24:04 2012 +0000 +++ b/main/queueing/inst/dtmc_fpt.m Fri Mar 16 14:08:44 2012 +0000 @@ -84,6 +84,10 @@ ( N>0 ) || \ error(err); + if ( any(diag(P) == 1) ) + error("Cannot compute first passage times for absorbing chains"); + endif + if ( nargin == 1 ) M = zeros(N,N); ## M(i,j) = 1 + sum_{k \neq j} P(i,k) M(k,j) @@ -119,6 +123,13 @@ %! 0.1 0.0 0.9; \ %! 0.9 0.1 0.0 ]; %! M = dtmc_fpt(P); +%! w = dtmc(P); +%! N = rows(P); +%! W = repmat(w,N,1); +%! Z = inv(eye(N)-P+W); +%! M1 = (repmat(diag(Z)',1,N) - Z) ./ repmat(w',1,N); +%! assert(M, M1); + %!shared P %! P = [ 0.0 0.9 0.1; \ @@ -138,10 +149,9 @@ %!test %! m = dtmc_fpt(P, 1, [2 3]); -## FIXME: check this (matrix not ergodic???) -%!xtest +%!test %! P = dtmc_bd([1 1 1], [ 0 0 0] ); -%! dtmc_fpt(P); +%! fail( "dtmc_fpt(P)", "absorbing" ); ## Example on p. 461 of ## http://www.cs.virginia.edu/~gfx/Courses/2006/DataDriven/bib/texsyn/Chapter11.pdf
--- a/main/queueing/inst/dtmc_is_irreducible.m Fri Mar 16 11:24:04 2012 +0000 +++ b/main/queueing/inst/dtmc_is_irreducible.m Fri Mar 16 14:08:44 2012 +0000 @@ -21,10 +21,11 @@ ## ## @cindex Markov chain, discrete time ## @cindex Discrete time Markov chain +## @cindex Irreducible Markov chain ## -## Check if @var{P} is irreducible, and identify strongly connected -## components in the transition graph of the discrete-time Markov chain -## with transition probability matrix @var{P}. +## Check if @var{P} is irreducible, and identify Strongly Connected +## Components (SCC) in the transition graph of the DTMC with transition +## probability matrix @var{P}. ## ## @strong{INPUTS} ## @@ -45,10 +46,10 @@ ## 1 if @var{P} irreducible, 0 otherwise. ## ## @item s -## @code{@var{s}(i)} is the strongly connected component that state @math{i} -## belongs to. Strongly connected components are numbered as 1, 2, -## @dots{}. If the graph is strongly connected, then there is a single -## SCC and the predicate @code{all(s == 1)} evaluates to true. +## @code{@var{s}(i)} is the SCC that state @math{i} belongs to. SCCs are +## numbered as 1, 2, @dots{}. If the graph is strongly connected, then +## there is a single SCC and the predicate @code{all(s == 1)} evaluates +## to true. ## ## @end table ## @@ -59,19 +60,23 @@ function [r s] = dtmc_is_irreducible( P ) - persistent epsilon = 10*eps; - if ( nargin != 1 ) print_usage(); endif - # dtmc_check_P(P); - + [N err] = dtmc_check_P(P); + if ( N == 0 ) + error(err); + endif s = __scc(P); r = (max(s) == 1); endfunction %!test +%! P = [0 .5 0; 0 0 0]; +%! fail( "dtmc_is_irresudible(P)" ); + +%!test %! P = [0 1 0; 0 .5 .5; 0 1 0]; %! [r s] = dtmc_is_irreducible(P); %! assert( r == 0 );
--- a/main/queueing/inst/dtmc_mtta.m Fri Mar 16 11:24:04 2012 +0000 +++ b/main/queueing/inst/dtmc_mtta.m Fri Mar 16 14:08:44 2012 +0000 @@ -22,6 +22,7 @@ ## ## @cindex Markov chain, disctete time ## @cindex Mean time to absorption +## @cindex Absorption probabilities ## ## Compute the expected number of steps before absorption for the ## DTMC with @math{N \times N} transition probability matrix @var{P}. @@ -53,9 +54,10 @@ ## @item B ## When called with a single argument, @var{B} is a @math{N \times N} ## matrix where @code{@var{B}(i,j)} is the probability of being absorbed -## in state @math{j}, starting from state @math{i}; if @math{j} is not -## absorbing, @code{@var{B}(i,j) = 0}; if @math{i} is absorbing, then -## @code{@var{B}(i,i) = 1}. When called with two arguments, @var{B} is a +## in state @math{j}, starting from transient state @math{i}; if +## @math{j} is not absorbing, @code{@var{B}(i,j) = 0}; if @math{i} is +## absorbing, then @code{@var{B}(i,i) = 1}. +## When called with two arguments, @var{B} is a ## vector with @math{N} elements where @code{@var{B}(j)} is the ## probability of being absorbed in state @var{j}, given initial state ## occupancy probabilities @var{p0}. @@ -82,7 +84,7 @@ if ( nargin == 2 ) ( isvector(p0) && length(p0) == K && all(p0>=0) && abs(sum(p0)-1.0)<epsilon ) || \ - usage( "p0 must be a probability vector" ); + usage( "p0 must be a state occupancy probability vector" ); endif