| Healthy | CVD | Death_CVD | Death_Other | |
|---|---|---|---|---|
| Healthy | 0 | 0 | 0 | 0 |
| CVD | 0 | 0 | 0 | 0 |
| Death_CVD | 0 | 0 | 0 | 0 |
| Death_Other | 0 | 0 | 0 | 0 |
Structuring and Adapting Models
iHEA Congress, July 2023
Motivation
Goal
A Novel Decision Model for Cardiovascular Disease (CVD) Prevention and Treatment in South Africa
Goal
An Adapted Decision Model for Cardiovascular Disease (CVD) Prevention and Treatment in South Africa
What Never Happens
How can we deal with this?
Solution
- Need “common ground” where we can combine and transform different model inputs.
- This “common ground” is often found in a transition rate matrix.
Solution: New Model Creation
Solution: New Model Creation
Solution: New Model Creation
Solution: Model Adaptation
Transition Rate Matrix
- The central “hub” of a Markov model.
- Straightforward to convert rate matrix into a transition probability matrix.
- Can be used to add transition states to record years of life lost (YLL) to disease.
- Facilitates direct modeling of DALYs–even in a model originally structured for QALYs!
Transition Rate Matrix
- Can also be used to change the cycle length.
- Facilitates modeling using alternative techniques:
- Continuous time Markov
- Discrete event simulation
- Microsimulation
Learning Objectives
- Understand differences between rates and probabilities, and rate-to-probability conversion formulas.
- Explain how to embed a transition probability matrix with a defined timestep.
- Demonstrate approaches for backwards conversion of an existing model.
I. Rates and Probabilities
What is the difference between a rate and a probability?
Rates vs. Probabilities
- Rate: number of occurrences of an event per unit of time.
- Probability: Liklihood that an event will occur for an in individual over a defined time period.
- Major difference is in denominator: rates take into account time at risk while probabilities do not.
Rates vs. Probabilities
Rates
- Number of events divided by the total time at risk experienced by all people followed.
- Ranges from 0 to ∞.
- #events in time periodTotal time period experienced by all subjects followed
Rates
Probabilities
Rates vs. Probabilities: Example1
- Suppose a study followed 100 people with cardiovascular disease for 4 years.
- At the end of 4 years, 40 had died.
- The probability of death over 4 years is 40/100=0.40.
Rates vs. Probabilities: Example1
- A rate takes into acccount the time each person was at risk.
- The 60 who survived were at risk the entire 4 years and contributed 60×4=240 years at risk.
- Once a person dies, he/she is no longer at risk.
Rates vs. Probabilities: Example1
- When a study does not report time at risk, the conventional assumption is that the events were spread evenly over the time period.
- Using this assumption, the average time at risk for the 40 who died was 2 years, adding 40 × 2 = 80 years at risk.
Rates vs. Probabilities: Example1
- Total time at risk for the cohort of 100 people is 320 person-years (240+80).
- Thus, the rate of death from CVD is 40/320 = 0.125 deaths per person-year.
- Pin this rate (0.125) into your mind for a few minutes; we’ll come back to it!
Summary: Rates vs. Probabilities
| Statistic | Evaluates | Range | Applicable Domain |
|---|---|---|---|
| Rate | #events in time periodTotal time period experienced by all subjects followed | 0 to ∞ | Rate matrix |
| Probability/risk | #events in time period#people followed for time period | 0-1 | Probability matrix |
II. Embedding a Model
What You Were Taught
- You were likely taught to convert each rate to a probability, then insert that probability into a transition matrix.
- This approach will yield the correct transition probabilities only if there is a single, non competing state you can transition to.
What We’ll Teach You
- We’ll teach you how to “embed” a transition probability matrix from a transition rate matrix.
- Essentially, fill out the transition matrix in terms of rates, then do a simple math operation to the entire matrix.
- This process is the matrix analogue to the usual rate-to-probability conversion formula.
Why Is This Important?
From the rate matrix we can…
- Swap in country-specific disease incidence and death rates, non-cause (“background”) mortality rates, etc.
- Add non-markovian transition states to calculate years of life lost to disease (YLL), a key input into DALYs.
- Add additional health states or strategies based on new evidence or country-specific context.
- Change the model time step (e.g., convert from monthly to yearly cycles).
Embedding 101
To construct a new model OR adapt an existing one, we must work from a transition rate matrix.
- “Embedding” refers to the process of converting a transition rate matrix to a transition probability matrix.
- We embed a transition probability matrix to using a four-step process.
Embedding: A Four-Step Process
- Use data inputs/published literature to define a rate matrix.
- Make strategy-specific adjustments to the rate matrix as needed.
- “Embed” the transition probability matrix using the rate matrix.
- Make further overall or strategy-specific adjustments to the transition probability matrix as needed.
1. Place rates in a rate matrix.
1. Place rates in a rate matrix.
- Let’s build on the earlier example to construct a rate matrix for a model with four health states: Healthy, CVD, Death from CVD, and Death from other causes.
1. Place rates in a rate matrix.
1. Place rates in a rate matrix.
- Recall that earlier, we calculated that the rate of death from CVD is 40/320 = 0.125 deaths per person-year.
- 0.125 per person year is the rate we’d enter into our rate matrix.
- Let’s also assume (for now) that the annual rate of CVD incidence in the population is 0.05.
1. Place rates in a rate matrix.
- Recall that earlier, we calculated that the rate of death from CVD is 40/320 = 0.125 deaths per person-year.
- 0.125 per person year is the rate we’d enter into our rate matrix.
- Let’s also assume (for now) that the annual rate of CVD incidence in the population is 0.05.
- We could also think of a separate death rate from background causes (e.g., 0.006 per person-year).1
1. Place rates in a rate matrix.
| Healthy | CVD | Death_CVD | Death_Other | |
|---|---|---|---|---|
| Healthy | 0 | 0.05 | 0.000 | 0.006 |
| CVD | 0 | 0.00 | 0.125 | 0.006 |
| Death_CVD | 0 | 0.00 | 0.000 | 0.000 |
| Death_Other | 0 | 0.00 | 0.000 | 0.000 |
1. Place rates in a rate matrix.
- The diagonal elements in a rate matrix are just the negative sum of the cells in the same row.
- In other words, the rows of a rate matrix sum to zero.
- Contrast this with a transition probability matrix, where the rows sum to one.
1. Place rates in a rate matrix.
- The diagonal elements in a rate matrix are just the negative sum of the cells in the same row.
- In this example, the diagonal value for the second row would be -0.131 = -(0.125+0.006)
- We can leave all other rate matrix values at 0 because the rate of progression from death to other states is zero.
- We’ve just completed step 1 of the four-step embedding process!
| Healthy | CVD | Death_CVD | Death_Other | |
|---|---|---|---|---|
| Healthy | -0.056 | 0.05 | 0.000 | 0.006 |
| CVD | 0 | -0.131 | 0.125 | 0.006 |
| Death_CVD | 0 | 0 | 0.000 | 0.000 |
| Death_Other | 0 | 0 | 0.000 | 0.000 |
1. Place rates in a rate matrix.
1. Place rates in a rate matrix.
- The rate matrix we just constructed represents a “natural history” model of the disease (CVD).
- This is a version of the model in which we allow the disease process to play out naturally, with no further intervention.
1. Place rates in a rate matrix.
- The “natural history” model is useful because it can help us verify that the model matches what we see in the real world.
- The “natural history” model also can be used to adapt transition rates, strategies, and other parameters to different countries/contexts.
For example, we would use the rate matrix to “swap in” a different background mortality rate, or a different disease incidence rate from another country, region, or population of interest.
1. Place rates in a rate matrix
- We’ll do this in the case study this afternoon when we draw on Global Burden of Disease (GBD) incidence and mortality data to “swap” in CVD rates for any country in the GBD data.
2. Make adjustments to the rate matrix as needed.
2. Make adjustments to the rate matrix as needed.
Suppose we want to:
- Adapt an existing model to a new country.
- Add strategies based on new evidence, or adapt strategies or parameters based on a decision problem in an LMIC.
2. Make adjustments to the rate matrix as needed.
- Think of the rate matrix as an automobile body shop: it’s where you go to make structural changes to your model.
- For now, we’ll assume you’re building a model from scratch.
- We’ll cover how to backwards-convert an existing model (based on a transition probability matrix) to a rate matrix later in this lecture.
2. Make adjustments to the rate matrix as needed.
- Suppose we want to model different strategies under consideration.
- We’ll do this by adapting a new rate matrix for each strategy.
- Note that not all adaptions can occur on the rate matrix; some may require adaptation of the transition probability matrix instead.
2. Make adjustments to the rate matrix as needed.
2. Make adjustments to the rate matrix as needed.
- Suppose that a new treatment (“Strategy A”) reduces the risk of CVD mortality by 20% (i.e., hazard ratio = 0.8).
- We can simply apply this hazard ratio directly to construct a rate matrix for strategy A.
- For the Strategy A rate matrix, the rate of CVD death is 0.8∗0.125=0.1
- Make sure that the diagonal element is adjusted to account for this change!
“Natural History” rate matrix
| Healthy | CVD | Death_CVD | Death_Other | |
|---|---|---|---|---|
| Healthy | -0.056 | 0.050 | 0.000 | 0.006 |
| CVD | 0.000 | -0.131 | 0.125 | 0.006 |
| Death_CVD | 0.000 | 0.000 | 0.000 | 0.000 |
| Death_Other | 0.000 | 0.000 | 0.000 | 0.000 |
2. Make adjustments to the rate matrix as needed.
- Suppose that a new treatment (“Strategy A”) reduces the risk of CVD mortality by 20% (i.e., hazard ratio = 0.8).
- We can simply apply this hazard ratio directly to construct a rate matrix for strategy A.
- For the Strategy A rate matrix, the rate of CVD death is 0.8∗0.125=0.1
- Make sure that the diagonal element is adjusted to account for this change!
“Natural History” rate matrix
| Healthy | CVD | Death_CVD | Death_Other | |
|---|---|---|---|---|
| Healthy | -0.056 | 0.050 | 0.000 | 0.006 |
| CVD | 0.000 | -0.131 | 0.125 | 0.006 |
| Death_CVD | 0.000 | 0.000 | 0.000 | 0.000 |
| Death_Other | 0.000 | 0.000 | 0.000 | 0.000 |
“Strategy A” rate matrix
| Healthy | CVD | Death_CVD | Death_Other | |
|---|---|---|---|---|
| Healthy | -0.056 | 0.05 | 0 | 0.006 |
| CVD | 0.000 | -0.106 | 0.1 | 0.006 |
| Death_CVD | 0.000 | 0 | 0 | 0.000 |
| Death_Other | 0.000 | 0 | 0 | 0.000 |
2. Make adjustments to the rate matrix as needed.
- Later today, we’ll show you how to further adapt the rate matrix to facilitate direct calculation of years of life lost to disease (YLL), a key input into DALY calculations.
3. “Embed” the transition probability matrix using the rate matrix
3. “Embed” the transition probability matrix using the rate matrix
3. “Embed” the transition probability matrix using the rate matrix
- Our next step is to convert the transition rate matrix into a transition probability matrix.
- Common practice is to use rate-to-probability conversion formulas:
3. “Embed” the transition probability matrix using the rate matrix
- Our next step is to convert the transition rate matrix into a transition probability matrix.
- Common practice is to use rate-to-probability conversion formulas:
3. “Embed” the transition probability matrix using the rate matrix
- Our next step is to convert the transition rate matrix into a transition probability matrix.
- Common practice is to use rate-to-probability conversion formulas
p=1−exp(−rt) where r is the rate and t is the time-step.
3. “Embed” the transition probability matrix using the rate matrix
- This formula works fine when there is only one possible state an individual can transition to.
- The formula will not calculate the correct transition probability if there are two or more states someone can transition to.
3. “Embed” the transition probability matrix using the rate matrix
- The technically correct approach for “embedding” a transition probability matrix is using the rate matrix exponential (Iosifescu 1980; Caswell and van Daalen 2021).
3. “Embed” the transition probability matrix using the rate matrix
- The technically correct approach for “embedding” a transition probability matrix is using the rate matrix exponential (Iosifescu 1980; Caswell and van Daalen 2021).
P=eRt
3. “Embed” the transition probability matrix using the rate matrix
- The technically correct approach for “embedding” a transition probability matrix is using the rate matrix exponential (Iosifescu 1980; Caswell and van Daalen 2021).
- This is the matrix analogue to the rate-to-probability conversion formula.
P=eRt
3. Calculate the Exponential of the Rate Matrix
P=eRt
3. Calculate the Exponential of the Rate Matrix
P=eRt
Pros:
Embedding the transition probability matrix in this way ensures that the correct transition probabilities are calculated.
This ensures that our discrete time Markov model accurately represents the underlying continuous time process.
3. Calculate the Exponential of the Rate Matrix
Cons:
- The computation of the matrix exponential requires complex numerical analysis methods (Moler and Van Loan 2003).
- Modern statistical software can now do it easily.
- We’ll show you a “hack” to do it in Excel!
3. Calculate the Exponential of the Rate Matrix
Cons:
- The computation of the matrix exponential requires complex numerical analysis methods (Moler and Van Loan 2003).
- Modern statistical software can now do it easily.
- We’ll show you a “hack” to do it in Excel!
3. Calculate the Exponential of the Rate Matrix
For our CVD example, the matrix exponential of the transition rate matrix yields the correct transition probability matrix:
| Healthy | CVD | Death_CVD | Death_Other | |
|---|---|---|---|---|
| Healthy | 0.9455 | 0.0455 | 0.0029 | 0.0060 |
| CVD | 0.0000 | 0.8772 | 0.1172 | 0.0056 |
| Death_CVD | 0.0000 | 0.0000 | 1.0000 | 0.0000 |
| Death_Other | 0.0000 | 0.0000 | 0.0000 | 1.0000 |
3. “Embed” the transition probability matrix using the rate matrix
- This process takes care of “combound” transitions and will result in some transition probabilities calculated that are seemingly impossible:
- Example: a non-zero probability of transition directly from a “Healthy” state to death from CVD in a single cycle.
Why do we obtain this probability?
| Healthy | CVD | Death_CVD | Death_Other | |
|---|---|---|---|---|
| Healthy | 0.946 | 0.046 | 0.003 | 0.006 |
| CVD | 0.000 | 0.877 | 0.117 | 0.006 |
| Death_CVD | 0.000 | 0.000 | 1.000 | 0.000 |
| Death_Other | 0.000 | 0.000 | 0.000 | 1.000 |
So Where Are We Now?
So Where Are We Now?
By constructing our model using the “roots” of a transition rate matrix, we can incorporate disparate sources of information.
Facilitates country/region-specific background mortality rates.
- Just swap in South Africa’s background mortality rate into the rate matrix, then embed the transition probability matrix!
- Or Columbia’s, or Thailand’s, …
- You could also swap in different disease incidence rates …
A Note of Caution
Not all literature-based parameters operate on transition rates.
You may find that the strategies you want to model have inputs based on odds ratios, relative risks, risk differences, etc.
- In these cases, you may need to embed the transition probability matrix first, then construct further adaptations from there.
4. Make adjustments to P as needed.
4. Make adjustments to P as needed.
- We must be careful about the scale on which model parameters apply.
- Odds ratios, relative risks, and risk differences all operate on the probability scale.
- That means you may need to make further adaptations to the transition probability matrix after you embed it.
4. Make adjustments to P as needed.
| Statistic | Evaluates | Range | Applicable Domain |
|---|---|---|---|
| Rate | #events in time periodTotal time period experienced by all subjects followed | 0 to ∞ | Rate matrix |
| Hazard Ratio | Hazard rate of outcome in exposedHazard rate of outcome in unexposed | 0 to ∞ | Rate matrix |
| Probability/risk | #events in time period#people followed for time period | 0-1 | Probability matrix |
| Odds | Probability of Outcome1−Probability of Outcome | 0 to ∞ | Probability matrix |
| Odds Ratio | Odds of outcome in exposedOdds of outcome in unexposed | 0 to ∞ | Probability matrix |
| Relative Risk | Probability of outcome in exposedProbablity of outcome in unexposed | 0 to ∞ | Probability matrix |
| Risk Difference | Probability of outcome in exposed−Probablity of outcome in unexposed | -1 to 1 | Probability matrix |
Full Process
III. Backwards Conversion
Goal
An Adapted Decision Model for HIV Prevention and Treatment in South Africa
III. Backwards Conversion
- What if we have a model that is already defined in terms of a transition probability matrix?
- How can we convert to a new country or setting?
- As shown earlier, we need the transition rate matrix to swap in new background mortality, new strategies, country-specific incidence, etc.
III. Backwards Conversion
III. Backwards Conversion
- We will next show you several ways to work backwards.
- Boils down to solving for the continuous “generator matrix” for the observed transition probability matrix.
- We’ll show you one approach to do this, which can be done in Excel but is technically incorrect (it will hopefully get you close, however!)
Working Example: HIV Model
Transition probability matrix:
| Healthy | LowCD4 | AIDS | Death | |
|---|---|---|---|---|
| Healthy | 0.721 | 0.202 | 0.067 | 0.010 |
| LowCD4 | 0.000 | 0.581 | 0.407 | 0.012 |
| AIDS | 0.000 | 0.000 | 0.750 | 0.250 |
| Death | 0.000 | 0.000 | 0.000 | 1.000 |
- From ‘Decision Modelling for Health Economic Evaluation’ by Briggs, Claxton, Sculpher.
Three Methods to Solve for the Generator
A. Convert the supplied transition probabilities to rates one-by-one (incorrect, but will get you close if you’re lucky!).
B. Compute the matrix logarithm of the transition probability matrix P.
C. Maximum likelihood-based approach.
You can always try 2+ of these methods to see how they compare.
A. Cellwise Conversion Using Rate-to-Probability Formula
- Simply convert each transition probability to a rate using standard rate-to-probability formula.
A. Cellwise Conversion Using Rate-to-Probability Formula
- Simply convert each transition probability to a rate using standard rate-to-probability formula.
r=−ln(1−p)
A. Cellwise Conversion Using Rate-to-Probability Formula
- Simply convert each transition probability to a rate using standard rate-to-probability formula.
- Note: technically incorrect, for all the reasons covered earlier relating to competing events and “jumpover” (i.e., compound transitions).
- Ideally, the underlying model would supply some underlying rates and you’d only have to do this for select transitions.
r=−ln(1−p)
A. Cellwise Conversion Using Rate-to-Probability Formula
Converted rate matrix:
| Healthy | LowCD4 | AIDS | Death | |
|---|---|---|---|---|
| Healthy | -0.305 | 0.226 | 0.069 | 0.010 |
| LowCD4 | 0.000 | -0.535 | 0.523 | 0.012 |
| AIDS | 0.000 | 0.000 | -0.288 | 0.288 |
| Death | 0.000 | 0.000 | 0.000 | 0.000 |
Converting our new rate matrix back to the probablity scale …
To see if this method gets us a rate matrix that is close to the generator matrix, we can use the methods discussed earlier to see if the rate matrix can recapitulate the original transition probability matrix…
Converting our new rate matrix back to the probablity scale …
Converted:
| Healthy | LowCD4 | AIDS | Death | |
|---|---|---|---|---|
| Healthy | 0.737 | 0.149 | 0.092 | 0.022 |
| LowCD4 | 0.000 | 0.586 | 0.347 | 0.067 |
| AIDS | 0.000 | 0.000 | 0.750 | 0.250 |
| Death | 0.000 | 0.000 | 0.000 | 1.000 |
Original:
| Healthy | LowCD4 | AIDS | Death | |
|---|---|---|---|---|
| Healthy | 0.721 | 0.202 | 0.067 | 0.010 |
| LowCD4 | 0.000 | 0.581 | 0.407 | 0.012 |
| AIDS | 0.000 | 0.000 | 0.750 | 0.250 |
| Death | 0.000 | 0.000 | 0.000 | 1.000 |
Converting our new rate matrix back to the probablity scale …
Converted:
| Healthy | LowCD4 | AIDS | Death | |
|---|---|---|---|---|
| Healthy | 0.737 | 0.149 | 0.092 | 0.022 |
| LowCD4 | 0.000 | 0.586 | 0.347 | 0.067 |
| AIDS | 0.000 | 0.000 | 0.750 | 0.250 |
| Death | 0.000 | 0.000 | 0.000 | 1.000 |
Original:
| Healthy | LowCD4 | AIDS | Death | |
|---|---|---|---|---|
| Healthy | 0.721 | 0.202 | 0.067 | 0.010 |
| LowCD4 | 0.000 | 0.581 | 0.407 | 0.012 |
| AIDS | 0.000 | 0.000 | 0.750 | 0.250 |
| Death | 0.000 | 0.000 | 0.000 | 1.000 |
Not great …
End of Main Slides (Optional Advanced Material Follows)
Advanced Approaches For Backwards Conversion
Alternative (Advanced Approaches)
- The two additional approaches (B. matrix logarithm; and C. maximum likelihood approach) will get you closer to a more accurate rate matrix.
- These approaches are beyond the capabilities of Excel.
- You will need to import the transition probability matrix into a statistical software package (e.g., R) to do them.
- We provide the process and code here, but will not cover it in person.
A Word of Warning
- As you’ll see, the next two approaches are not always going to work perfectly. Or they may require you to impose some structural assumptions on your model to make them work.
- If the original transition probability matrix was defined incorrectly (e.g., no jumpover probabilities), the generator matrix may not exist.
- Identifiability is a more general issue, however.
B. Matrix Logarithm of P
- In mathematical terms, the generator matrix is the matrix logarithm of the transition probability matrix.
- A matrix has a logarithm if and only if if it is invertible.
B. Matrix Logarithm of P
- In mathematical terms, the generator matrix is the matrix logarithm of the transition probability matrix.
- A matrix has a logarithm if and only if if it is invertible.
R=logP
B. Matrix Logarithm of P
- The log can be found using spectral or eigenvalue decomposition.
- This is beyond the capabilities of Excel; what follows is code for the open-source R language.
- If V is a matrix where each column is an eigenvector of P, then
B. Matrix Logarithm of P
- The log can be found using spectral or eigenvalue decomposition.
- This is beyond the capabilities of Excel; what follows is code for the open-source R language.
- If V is a matrix where each column is an eigenvector of P, then
A′=V−1AV
B. Matrix Logarithm of P
- The log can be found using spectral or eigenvalue decomposition.
- This is beyond the capabilities of Excel; what follows is code for the open-source R language.
- If V is a matrix where each column is an eigenvector of P, then
A′=V−1AV
logP=V(logA′)V−1
B. Matrix Logarithm of P
# read the transition probability matrix as stored in an
# excel file
mP <- readxl::read_xlsx("path-excel-file.xls")
# eigenvalue decomposition
V <- eigen(mP)$vectors
iV <- solve(V)
Ap <- iV %*% mP %*% V
lAp <- diag(log(diag(Ap)), nrow(Ap), ncol(Ap))
R <- V %*% lAp %*% iV
R[abs(R) < 1e-6 ] <- 0
dimnames(R) = dimnames(mP)
round(R,3)B. Matrix Logarithm of P
V <- eigen(mP)$vectors
iV <- solve(V)
Ap <- iV %*% mP %*% V
lAp <- diag(log(diag(Ap)), nrow(Ap), ncol(Ap))
R <- V %*% lAp %*% iV
R[abs(R) < 1e-6 ] <- 0
dimnames(R) = dimnames(mP)
round(R,3) Healthy LowCD4 AIDS Death
Healthy -0.327 0.311 0.002 0.013
LowCD4 0.000 -0.543 0.615 -0.072
AIDS 0.000 0.000 -0.288 0.288
Death 0.000 0.000 0.000 0.000B. Matrix Logarithm of P
| Healthy | LowCD4 | AIDS | Death | |
|---|---|---|---|---|
| Healthy | -0.327 | 0.311 | 0.002 | 0.013 |
| LowCD4 | 0 | -0.543 | 0.615 | -0.072 |
| AIDS | 0 | 0 | -0.288 | 0.288 |
| Death | 0 | 0 | 0 | 0 |
- Note there is a negative rate from LowCD4 Death!
- This is a transition from Death LowCD4 and should be moved to the other side of the matrix.
B. Matrix Logarithm of P
| Healthy | LowCD4 | AIDS | Death | |
|---|---|---|---|---|
| Healthy | -0.327 | 0.311 | 0.002 | 0.013 |
| LowCD4 | 0 | -0.543 | 0.615 | -0.072 |
| AIDS | 0 | 0 | -0.288 | 0.288 |
| Death | 0 | 0 | 0 | 0 |
- The negative rate implies an identifiability issue.
- Essentially, the math only works out if we allow for negative rates—in this case, bringing people back alive from the death state!
B. Matrix Logarithm of P
expm::logmHigham (2008) method returns same as eigenvalue method.- Tweaking rates to get a proper model is ad-hoc and difficult.
B. Matrix Logarithm of P
- Higham (2008) method returns same as eigenvalue method. -
expm::logm() - Tweaking rates to get a proper model is ad-hoc and difficult.
Healthy LowCD4 AIDS Death
Healthy -0.327 0.311 0.002 0.013
LowCD4 0.000 -0.543 0.615 -0.072
AIDS 0.000 0.000 -0.288 0.288
Death 0.000 0.000 0.000 0.000C. Multi-State Model Based Approach
- Imposition of structural assumptions and fitting to pseudo-data derived from original Markov model can result in reasonable rate estimates.
- We could assume that a patient with HIV at this point in time only gets sicker and external causes of death are negligible. Healthy LowCD4 AIDS Death constrains model.
- We use the reported data from the original model at year 1.
C. Multi-State Model Based Approach
Steps:
- Run the original model for 1+ cycles to obtain the Markov trace.
- Construct psuedo-data from the resulting trace based on a cohort with reasonable size (e.g., 1,000 patients)
- Estimate transition hazards in the pseudo-data based on a multistate model.1
- Use the estimated transition hazards to construct the rate matrix.
1. Run the original model
tr0 <- # Initial state occupancy.
c("Healthy" = 1000,"LowCD4" = 0, "AIDS" = 0, "Death" = 0)
tr1 <- # State occupancy after one cycle.
tr0 %*% mP
tr <- # Bind together into a data frame.
rbind.data.frame(tr0, tr1) %>%
mutate(cycle = c(0, 1)) %>%
select(cycle,everything())
tr cycle Healthy LowCD4 AIDS Death
1 0 1000 0 0 0
2 1 721 202 67 102. Construct psuedo-data
- Idea: create data for 1,000 simulated “patients” followed for two cycles.
- Counts in the Markov trace govern state occupancy in each cycle.
- In first cycle, all 1,000 are in
Healthystate. - In second cycle, 721 remain in
Healthywhile 202 are inLowCD4, etc.
- In first cycle, all 1,000 are in
2. Construct psuedo-data
tr <-
rbind.data.frame(tr0, tr1) %>%
mutate(cycle = c(0, 1)) %>%
gather(state,count,-cycle) %>%
mutate(count = round(count,0)) %>%
arrange(cycle) %>%
lapply(.,rep,.$count) %>%
cbind.data.frame() %>%
as_tibble() %>%
mutate(state = factor(state,
levels = c("Healthy","LowCD4","AIDS","Death"))) %>%
arrange(cycle,state) %>%
mutate(ptnum = rep(1:1000,2)) %>%
select(-count) %>%
arrange(ptnum,cycle) %>%
mutate(state = as.numeric(state)) %>%
select(ptnum,cycle,state)# A tibble: 2,000 × 3
ptnum cycle state
<int> <dbl> <dbl>
1 1 0 1
2 1 1 1
3 2 0 1
4 2 1 1
5 3 0 1
6 3 1 1
7 4 0 1
8 4 1 1
9 5 0 1
10 5 1 1
# ℹ 1,990 more rows3. Fit Multistate Model
- Next, fit a multistate model to the pseudo-data.
- The multistate model will use likelihood-based methods to estimate transition intensities among health states in the data.
- These transition intensities are the rates that you can use in the rate matrix!
3. Fit Multistate Model
- Define the state table for the multistate model
- Check that cell transition counts match the trace.
3. Fit Multistate Model
- Initial rate guesses are based on the eigenvalue method, i.e., result of
logm(mP).
- Fit the model
3. Fit Multistate Model
Call:
msm(formula = state ~ cycle, subject = ptnum, data = tr, qmatrix = Q.init)
Maximum likelihood estimates
Transition intensities
Baseline
Healthy - Healthy -0.32712 ( -3.680e-01, -2.907e-01)
Healthy - LowCD4 0.32702 ( 1.354e-01, 7.898e-01)
Healthy - Death 0.00010 ( 0.000e+00, Inf)
LowCD4 - LowCD4 -0.64477 ( -1.139e+01, -3.650e-02)
LowCD4 - AIDS 0.63465 ( 9.115e-06, 4.419e+04)
LowCD4 - Death 0.01012 ( 8.238e-299, 1.242e+294)
AIDS - AIDS -0.34837 ( -3.805e+40, -3.189e-42)
AIDS - Death 0.34837 ( 3.189e-42, 3.805e+40)
-2 * log-likelihood: 1572.2 4. Construct the Rate Matrix
5. Embed the Transition Probability Matrix
Multistate-Model:
Markov Trace Comparison
Healthy LowCD4 AIDS Death
[1,] 721 202 67 10 Healthy LowCD4 AIDS Death
[1,] 721 202 67 10- Resulting expectation after 1 cycle is nearly identical.
- The msm method has bits down at 4th decimal. Round off made them identical.
- Only a tiny bit of error (721 vs. 721.0005) makes
lognumerically unstable converting from probability to rate.
Key Takeaways
- By imposing some constraints on the underlying transitions, we were able to yield a generator matrix that makes sense!
- Critically, our multistate model-based generator matrix closely approximates the original transition probability matrix, but does not imply negative rates that bring people back from death.
Key Takeaways
Multistate-Model:
Key Takeaways
With a reasonable generator matrix defined, we can now augment the model as we see fit:
- Add additional health states.
- Add new strategies with evidence on efficiacy drawn from the literature (e.g., hazard rates).
- Add accumulators and transition states.
- Change the cycle length.
Once done with the above, just embed the new matrix.
References
Caswell, Hal, and Silke van Daalen. 2021. “Healthy Longevity from Incidence-Based Models.” Demographic Research 45: 397–452.
Gidwani, Risha, and Louise B. Russell. 2020. “Estimating Transition Probabilities from Published Evidence: A Tutorial for Decision Modelers.” PharmacoEconomics 38 (11): 1153–64. https://doi.org/10.1007/s40273-020-00937-z.
Higham, Nicholas J. 2008. Functions of Matrices: Theory and Computation. SIAM.
Iosifescu, M. 1980. “Finite Markov Processes and Their Applications. Wiley.” New York.
Moler, Cleve, and Charles Van Loan. 2003. “Nineteen Dubious Ways to Compute the Exponential of a Matrix, Twenty-Five Years Later.” SIAM Review 45 (1): 3–49. https://doi.org/10.1137/S00361445024180.
Structuring and Adapting Models iHEA Congress, July 2023