# Markov Chain Transition Intensity Conversion

I have a question about converting a 3-state discrete state, continuous-time, markov chain to a 2-state.

My 3-state model has states: Well (state 1), Ill (state 2) and Dead (state 3).

$$\begin{bmatrix}-(a12 + a13) & a12 & a13\\0 & -a23 & a23\\ 0 & 0 & 0\end{bmatrix}$$ This 3 state matrix is full.mat in the R code.

I would like to convert it to an Alive/Dead model. I am not sure if I can do it the following way: $$\begin{bmatrix}-(a13 + a23) & (a13+a23)\\0 & 0\end{bmatrix}$$ where I am simply adding the intensity of Well->Dead, and Ill->Dead to compute the intensity of Alive->Dead for a 2-state model? This matrix is small.mat in the R code.

I would expect that the sum of transition probabilities P(1->3) + P(2->3) from the three state model should equal P(alive -> dead) in the 2-state model.

Essentially, I am trying to determine $$\mathrm{Pr}(X(t+h) = 3 | X(t) =1~~OR~~X(t) =2)$$ But the final 2 lines of the R-code show that these values are not equivalent, they are slightly off... Am I doing things incorrectly, or is this just rounding approximation by expm()?

library(expm)

full.mat<- rbind(c(-0.003260632, 0.000514263, 0.002746369),
c(0.000000000, -0.007948859, 0.007948859),
c(0.000000000, 0.000000000, 0.000000000))

small.mat<-matrix(0,2,2)
small.mat[1,2]<-full.mat[1,3]+full.mat[2,3]
small.mat[1,1]<-small.mat[1,2]*-1

exp.full<-expm(full.mat)
exp.small<-expm(small.mat)

# COMPUTE PROBABILITY OF DEATH
exp.small[1,2] # this is probability of death in 2-state model
exp.full[1,3]+exp.full[2,3] # this is probability of death
in 3-state model

-

The first thing to realize is that in this model one needs to know whether one is well-being or ill to know the "chances" one has to become dead.

The exception is when $a_{13}=a_{23}=\alpha$, then the alive/dead process is indeed a Markov process on the state space $\{\mathtt{alive},\mathtt{dead}\}$ with rate transition matrix $\begin{pmatrix}-\alpha &\alpha\\ 0 & 0\end{pmatrix}$. In every other case, the usual Bayes decomposition yields $$\mathbb P(X(t+\mathrm dt)=3\mid X(t)=1\ \text{or}\ 2)=\alpha(t)\mathrm dt,$$ where $$\alpha(t)=\frac{a_{13}p_1(t)+a_{23}p_2(t)}{p_1(t)+p_2(t)},\qquad p_i(t)=\mathbb P(X(t)=i).$$ Note that each $p_i(t)$ depends on $t$ and on the initial distribution $(p_1(0),p_2(0))$. Recall that $$p_1(t)=p_1(0)\mathrm e^{-(a_{12}+a_{13})t},\quad p_2(t)=p_1(0)c(t)+p_2(0)\mathrm e^{-a_{23}t},$$ for some explicit function $c(t)$ you might want to write down.

To sum up, call $Y(t)=\mathtt{dead}$ if $X(t)=3$ and $Y(t)=\mathtt{alive}$ otherwise. Then $(Y(t))_{t\geqslant0}$ is not (in general) a Markov process on the state space $\{\mathtt{alive},\mathtt{dead}\}$ because the distribution of the state $Y(t+\mathrm dt)$ depends on the distribution of the state $Y(t)$ (good), on $t$ itself (medium good), and also on the initial distribution of $X(0)$ (not good).

-
Thank you. So because this isn't a markov process, I can't use the matrix exponential of the transition intensity matrix to compute the transition probability matrix? – MPahuta Feb 3 '13 at 21:34
There is no transition intensity matrix, if you think about it (except in the singular case $a_{13}=a_{23}$). – Did Feb 4 '13 at 12:11

All rates defining the progressive model can be organized in a rate matrix $\mathbf{\mathrm{R}}$. $\mathbf{\mathrm{R}}$ is organized such that rows correspond to the "from" state, and the columns to the "to" state. Disallowed transitions are assigned a rate of $0$, and each row must sum to $0$. \begin{align*} \mathbf{\mathrm{R}}=\left[\begin{array}{c c c} -(\alpha+b \alpha) & b \alpha &\alpha\\ 0 & -g \alpha & g \alpha\\ 0 & 0 & 0\\ \end{array} \right]. \end{align*} We also define the state occupation column vector $\mathbf{\mathrm{P}}(t)$ \begin{align*} \mathbf{\mathrm{P}}(t)=\left[\begin{array}{c} W(t)\\ I(t)\\ D_O(t)\\ \end{array} \right]. \end{align*} The rate matrix is used to compute the state occupation column vector at any arbitrary time $t$ using the relationship\cite{Cox1965} \begin{align} \mathbf{\mathrm{P}}(t) &= \mathbf{\mathrm{P}}(0) e^{\mathrm{\textbf{R}} t} \label{eq:p_rec}. \end{align} By adding two constraints, \begin{align} D_O(t) &= 1 - (W(t) + I(t)),~\mathrm{and}\\ D_O(0) &=0, \end{align} we can derive a general formula for the overall death function: $$D_O(t) = 1 - W(0) e^{- \alpha t (1+b)} - W(0) \frac{b }{g - 1- b }\left(e^{- \alpha t (1+b)} - e^{- \alpha t g}\right)- \left( 1 - W(0)\right)e^{-\alpha t g}$$

-