Time- and Country/Region-Varying Transition Probabilities

Introduction

Another problem that occurs is that of non-stationary rates. What if a rate is changing over an interval time step? How is this best modeled? A first order approach would be to simply take the mean of the rates at either time point. A higher order approach is to integrate the total risk exposure over the time interval.

\[r_{e} = \int_{t_1}^{t_2} r(t) dt\] This computes the effective rate \(r_e\) from a time varying rate, \(r(t)\).

For example, the probability of death from outside causes could be modeled based on mortality tables using a Gompertz distribution with a shape \(a\) and rate \(b\) parameter. The hazard function of a distribution is the continuous rate function. The hazard function can be derived from the probability density function (PDF) divided by 1 minus cumulative probability function (CDF). For the Gompertz this is,

\[h(t) = \frac{f(t)}{1-F(t)} = \frac{b e^{at} e^{-b/a(e^{at}-1)}}{e^{-b/a (e^{at}-1)}} = b e^{at}\]

From this the effective rate over a timestep of a discrete Markov simulation can be computed.

\[r_{e} = \int_{t_1}^{t_2} b e^{at} dt = \frac{b}{a} (e^{at_2} - e^{at_1})\]

While this exists as an anayltic result, in general it’s better to use the effective risk over a time interval as this leads to improved accuracy of a simulation.

Setup

Let’s start with the original rate matrix for the healthy-sick-dead model. Below is the rate matrix for the standard of care (“Treatment A”) strategy:

	Healthy	Sick	Dead
Healthy	-(r_HS + r_HD)	r_HS	r_HD
Sick	0	-(hr_S * r_HD)	hr_S * r_HD
Dead	0	0	0

Time- and Country/Region-Specific Death Rates

Now suppose that we can parameterize backround mortality using a shape parameter (\(a\)) and a scale parameter (\(b\)) estimated using life table data from our country/region of interest.

Age	shape	scale	cycle
25	0.08854	0.00050	0
26	0.08877	0.00054	1
27	0.08900	0.00058	2
28	0.08921	0.00063	3
29	0.08942	0.00068	4
30	0.08962	0.00074	5
75	0.09770	0.03911	50
76	0.09744	0.04336	51
77	0.09678	0.04831	52
78	0.09610	0.05381	53
79	0.09521	0.06001	54
80	0.09389	0.06723	55

We can calculate the specific mortality rate by age using the following function:

rate_secular_death_by_age <- function(t, h=1, a, b)
{
  (b/a) * (exp(a*t) - exp(a*(t-h)))
}

cycle	Age	shape	scale	r_HD_gompertz
0	25	0.08854	0.00050	0.00048
1	26	0.08877	0.00054	0.00051
2	27	0.08900	0.00058	0.00056
3	28	0.08921	0.00063	0.00060
4	29	0.08942	0.00068	0.00065
5	30	0.08962	0.00074	0.00070
50	75	0.09770	0.03911	0.03726
51	76	0.09744	0.04336	0.04131
52	77	0.09678	0.04831	0.04604
53	78	0.09610	0.05381	0.05130
54	79	0.09521	0.06001	0.05725
55	80	0.09389	0.06723	0.06417

Notice in the table that one column (r_HD_gompertz) provides the death rate for a given age. We’ll pull out this column and define it as a new vector parameter called r_HD_gompertz.

We next update our transition rate definitions to include this new parameter:

	Healthy	Sick	Dead
Healthy	-(r_HS + r_HD_gompertz[t])	r_HS	r_HD_gompertz[t]
Sick	0	-(hr_S * r_HD_gompertz[t])	hr_S * r_HD_gompertz[t]
Dead	0	0	0

We can now use this rate matrix in lieu of the fixed rate matrix we started with. Owing to time constraints in the workshop, we will not cover how to embed these matrices in Excel, however you can find a nice primer on how to create time-varying transition Markov models here.

We now have the time-varying transition probability matrices:

	Healthy	Sick	Dead
Healthy	0.8603	0.13916	0.00055
Sick	0.0000	0.99857	0.00143
Dead	0.0000	0.00000	1.00000

	Healthy	Sick	Dead
Healthy	0.43759	0.03818	0.52423
Sick	0.00000	0.13142	0.86858
Dead	0.00000	0.00000	1.00000

Let’s go one step further and intoduce age-dependence into the rate of progressing from sickness to death. Suppose a research study estimates the following age-stratified hazard ratios:

age	Hazard Ratio (Sickness to Death)
<30	1.5
30-49	2
50-74	3
75+	4

We can then define a parameter vector for the HR_S parameter (which adjusts the regular death rate upwards to account for the higher likelihood of death after becoming sick) based on the age in the model.

age	hr_S
25	1.5
26	1.5
27	1.5
28	1.5
29	1.5
30	2.0
74	3.0
75	4.0
76	4.0
77	4.0
78	4.0
79	4.0
80	4.0

Here is our amended rate matrix

	Healthy	Sick	Dead
Healthy	-(r_HS + r_HD_gompertz[t])	r_HS	r_HD_gompertz[t]
Sick	0	-(hr_S[t] * r_HD_gompertz[t])	hr_S[t] * r_HD_gompertz[t]
Dead	0	0	0

	Healthy	Sick	Dead
Healthy	0.8603	0.13921	0.00049
Sick	0.0000	0.99928	0.00072
Dead	0.0000	0.00000	1.00000

	Healthy	Sick	Dead
Healthy	0.43759	0.02959	0.53281
Sick	0.00000	0.06681	0.93319
Dead	0.00000	0.00000	1.00000

We now have a series of age-specific transition probability matrices for our decision problem. The rest of the analysis proceeds as normal, however we just need to make sure we are applying the correct transition probability matrix at each cycle.