aphp / heemod

Markov Models for Health Economic Evaluations
https://aphp.github.io/heemod/
GNU General Public License v3.0
11 stars 1 forks source link

dsa in heemod package #6

Closed huipharmacist closed 1 year ago

huipharmacist commented 1 year ago

When doing dsa, some parameters in the dsa diagram are not displayed. Why?thank you.

5

image

KZARCA commented 1 year ago

Hi, could you please give me a code example? Best

huipharmacist commented 1 year ago

define_dsa<- define_dsa(cLapa,26, 27.6,cCape, 1.2, 1.8,BSA, 1.5, 1.9, uS, 0.63, 0.77,uR, 0.7, 0.924,uP, 0.45, 0.55) result_dsa<- run_dsa(res_mod, dsa = define_dsa) plot(result_dsa, result = "icer", type = "difference")

The code should be fine. I want to ask that base ICER is not included between the minimum and maximum range of the cLapa parameter in the figure above, but it still extends to base ICER. So numbers are covered up. In cases like this, how to modify the code so that it shows the real results. I hope you can help me.

KZARCA commented 1 year ago

Hi, I need your res_mod object to help you. Can you provide me a simplified version of your code just to be able to reproduce the problem by myself? Best

huipharmacist commented 1 year ago

The previous example may not be so appropriate, it was written casually, please see my new one, the same problem, thank you

mod_par<-define_parameters(

pHY1toHY2=0.3604,

pHY2toHY2.5_AA=0.0654,

pHY2.5toHY3_AA=0.1758,

pHY3toHY4=0.2079,

pHY4toHY5=0.2092,

pHYtoDeath_AA=0.0347,

HRLDvsMAOBI=0.85,

pHY1toHY1_AA=1-pHY1toHY2-pHYtoDeath_AA,

pHY2toHY2_AA=1-pHY2toHY2.5_AA-pHYtoDeath_AA,

pHY2.5toHY2.5_AA=1-pHY2.5toHY3_AA-pHYtoDeath_AA,

pHY3toHY3_AA=1-pHY3toHY4-pHYtoDeath_AA,

pHY4toHY4_AA=1-pHY4toHY5-pHYtoDeath_AA,

pHY5toHY5_AA=1-pHYtoDeath_AA,

pHY2toHY2.5_BB=ifelse(markov_cycle<7,0.0654,0.0157),

pHY2.5toHY3_BB=ifelse(markov_cycle<7,0.1758,0.0615),

pHYtoDeath_BB= pHYtoDeath_AA/HRLDvsMAOBI,

pHY1toHY1_BB=1-pHY1toHY2-pHYtoDeath_BB,

pHY2toHY2_BB=1-pHY2toHY2.5_BB-pHYtoDeath_BB,

pHY2.5toHY2.5_BB=1-pHY2.5toHY3_BB-pHYtoDeath_BB,

pHY3toHY3_BB=1-pHY3toHY4-pHYtoDeath_BB,

pHY4toHY4_BB=1-pHY4toHY5-pHYtoDeath_BB,

pHY5toHY5_BB=1-pHYtoDeath_BB)

添加修改参数,健康效用值以及成本

mod_par<-modify(mod_par,

            cInspection=796.92,

            cHY1to2.5Medical=51.7,

            cHY3Medical=74.62,

            cHY4Medical=89.26,

            cHY5Medical=65.24,

            cHY1to2.5rehabilitation=33.03,

            cHY3rehabilitation=221.25,

            cHY4rehabilitation=441.02,

            cHY5rehabilitation=605.77,

            cComplications=785.57,

            cdrug_AA=1.93,

            cHY3motor_AA=8706.88,

            cHY4motor_AA=11229.02,

            cHY5motor_AA=11743.52,

            uHY1to2.5_AA=0.708360851355,

            uHY3to5_AA=0.658386842,

            cdrug_BB=3.18,

            cHY3motor_BB=6934.73,

            cHY4motor_BB=8905.51,

            cHY5motor_BB=9953.94,

            uHY1to2.5_BB=0.706802940955,

            uHY3to5_BB=0.65659185552,

            r=0.05)

mat_trans_AA<-define_transition(

state_names=c("HY1","HY2","HY2.5","HY3","HY4","HY5","Death"),

pHY1toHY1_AA,pHY1toHY2,0,0,0,0,pHYtoDeath_AA,

0,pHY2toHY2_AA,pHY2toHY2.5_AA,0,0,0,pHYtoDeath_AA,

0,0,pHY2.5toHY2.5_AA,pHY2.5toHY3_AA,0,0,pHYtoDeath_AA,

0,0,0,pHY3toHY3_AA,pHY3toHY4,0,pHYtoDeath_AA,

0,0,0,0,pHY4toHY4_AA,pHY4toHY5,pHYtoDeath_AA,

0,0,0,0,0,pHY5toHY5_AA,pHYtoDeath_AA,

0,0,0,0,0,0,1)

mat_trans_BB<-define_transition(

state_names=c("HY1","HY2","HY2.5","HY3","HY4","HY5","Death"),

pHY1toHY1_BB,pHY1toHY2,0,0,0,0,pHYtoDeath_BB,

0,pHY2toHY2_BB,pHY2toHY2.5_BB,0,0,0,pHYtoDeath_BB,

0,0,pHY2.5toHY2.5_BB,pHY2.5toHY3_BB,0,0,pHYtoDeath_BB,

0,0,0,pHY3toHY3_BB,pHY3toHY4,0,pHYtoDeath_BB,

0,0,0,0,pHY4toHY4_BB,pHY4toHY5,pHYtoDeath_BB,

0,0,0,0,0,pHY5toHY5_BB,pHYtoDeath_BB,

0,0,0,0,0,0,1)

plot(mat_trans_AA)

plot(mat_trans_BB)

par(mfrow = c(1,2))

3.定义各个状态下的参数值(主要是成本和健康效用值)

state_HY1<-define_state(

cost=discount(

dispatch_strategy(

  BB=ifelse(markov_cycle<1,cdrug_BB/5*(5*7+10*23*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications,

            cdrug_BB/5*10*24*7+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications),

  AA=ifelse(markov_cycle<1,cdrug_AA/250*((125+250+375+500+625)*7+750*19*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications,

            cdrug_AA/250*(750*24*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications)

),

r),

qaly=discount(

dispatch_strategy(

  AA=uHY1to2.5_AA,

  BB=uHY1to2.5_BB),

r)

)

state_HY2<-define_state(

cost=discount(

dispatch_strategy(

  AA=ifelse(markov_cycle<1,cdrug_AA/250*((125+250+375+500+625)*7+750*19*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications,

            cdrug_AA/250*(750*24*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications),

  BB=ifelse(markov_cycle<1,cdrug_BB/5*(5*7+10*23*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications,

            cdrug_BB/5*10*24*7+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications)),

r),

qaly=discount(

dispatch_strategy(

  AA=uHY1to2.5_AA,

  BB=uHY1to2.5_BB),

r)

)

state_HY2.5<-define_state(

cost=discount(

dispatch_strategy(

  AA=ifelse(markov_cycle<1,cdrug_AA/250*((125+250+375+500+625)*7+750*19*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications,

            cdrug_AA/250*(750*24*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications),

  BB=ifelse(markov_cycle<1,cdrug_BB/5*(5*7+10*23*7)+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications,

            cdrug_BB/5*10*24*7+cInspection+cHY1to2.5Medical+cHY1to2.5rehabilitation+cComplications)),

r),

qaly=discount(

dispatch_strategy(

  AA=uHY1to2.5_AA,

  BB=uHY1to2.5_BB),

r)

)

state_HY3<-define_state(

cost=discount(

dispatch_strategy(

  AA=ifelse(markov_cycle<1,cdrug_AA/250*((125+250+375+500+625)*7+750*19*7)+cHY3Medical+cHY3rehabilitation+cHY3motor_AA+cComplications,

            cdrug_AA/250*(750*24*7)+cHY3Medical+cHY3rehabilitation+cHY3motor_AA+cComplications),

  BB=ifelse(markov_cycle<1,cdrug_BB/5*(5*7+10*23*7)+cHY3Medical+cHY3rehabilitation+cHY3motor_BB+cComplications,

            cdrug_BB/5*10*24*7+cHY3Medical+cHY3rehabilitation+cHY3motor_BB+cComplications)),

r),

qaly=discount(

dispatch_strategy(

  AA=uHY3to5_AA,

  BB=uHY3to5_BB),

r)

)

state_HY4<-define_state(

cost=discount(

dispatch_strategy(

  AA=ifelse(markov_cycle<1,cdrug_AA/250*((125+250+375+500+625)*7+750*19*7)+cHY4Medical+cHY4rehabilitation+cHY4motor_AA+cComplications,

            cdrug_AA/250*(750*24*7)+cHY4Medical+cHY4rehabilitation+cHY4motor_AA+cComplications),

  BB=ifelse(markov_cycle<1,cdrug_BB/5*(5*7+10*23*7)+cHY4Medical+cHY4rehabilitation+cHY4motor_BB+cComplications,

            cdrug_BB/5*10*24*7+cHY4Medical+cHY4rehabilitation+cHY4motor_BB+cComplications)),

r),

qaly=discount(

dispatch_strategy(

  AA=uHY3to5_AA,

  BB=uHY3to5_BB),

r)

)

state_HY5<-define_state(

cost=discount(

dispatch_strategy(

  AA=ifelse(markov_cycle<1,cdrug_AA/250*((125+250+375+500+625)*7+750*19*7)+cHY5Medical+cHY5rehabilitation+cHY5motor_AA+cComplications,

            cdrug_AA/250*(750*24*7)+cHY5Medical+cHY5rehabilitation+cHY5motor_AA+cComplications),

  BB=ifelse(markov_cycle<1,cdrug_BB/5*(5*7+10*23*7)+cHY5Medical+cHY5rehabilitation+cHY5motor_BB+cComplications,

            cdrug_BB/5*10*24*7+cHY5Medical+cHY5rehabilitation+cHY5motor_BB+cComplications)),

r),

qaly=discount(

dispatch_strategy(

  AA=uHY3to5_AA,

  BB=uHY3to5_BB),

r)

)

state_Death<-define_state(

cost= 0,

qaly= 0)

strat_AA<-define_strategy(

transition=mat_trans_AA,

HY1=state_HY1,

HY2=state_HY2,

HY2.5=state_HY2.5,

HY3=state_HY3,

HY4=state_HY4,

HY5=state_HY5,

Death=state_Death)

strat_BB<-define_strategy(

transition=mat_trans_BB,

HY1=state_HY1,

HY2=state_HY2,

HY2.5=state_HY2.5,

HY3=state_HY3,

HY4=state_HY4,

HY5=state_HY5,

Death=state_Death)

res_mod<-run_model(

parameters=mod_par,

BB=strat_BB,

AA=strat_AA,

cycles= 100,

cost=cost,

effect=qaly,

method = "life-table")

res_mod

cycle_AA<-data.frame(res_mod$eval_strategy_list$AA$parameters)

cycle_BB<-data.frame(res_mod$eval_strategy_list$gBB$parameters)

plot(res_mod)

plot(res_mod,type="ce")

敏感性分析

单因素敏感性分析

def_dsa<-define_dsa(

pHY1toHY2,0.2523,0.4685,

pHY2toHY2.5_AA,0.0458,0.0850,

pHY2.5toHY3_AA,0.1231,0.2285,

pHY3toHY4,0.1455,0.2703,

pHY4toHY5,0.1464,0.2720,

HRLDvsMAOBI,0.69,1.06,

r,0,0.08,

uHY1to2.5_AA,0.6862,0.7283,

uHY3to5_AA,0.6329,0.6815,

uHY1to2.5_BB,0.6992,0.7030,

uHY3to5_BB,0.6502,0.6560

)

res_dsa<-run_dsa(res_mod,def_dsa)

options(max.print = 10000)

res_dsa

plot(res_dsa,result="icer",type="difference")

At 2023-06-01 19:59:30, "KZARCA" @.***> wrote:

Hi, I need your res_mod object to help you. Can you provide me a simplified version of your code just to be able to reproduce the problem by myself? Best

— Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you authored the thread.Message ID: @.***>

huipharmacist commented 1 year ago

Is there any way to solve it, thank you.

KZARCA commented 1 year ago

The problem I found is that the ICER for the 2 extreme values of your parameter does not contain the baseline ICER. I've never seen this before. It is possible that your code contains an error.

huipharmacist commented 1 year ago

thank you, is there any way to adjust the graph manually through ggplot

KZARCA commented 1 year ago

Hi, the function is plot.dsa; what you could do is creating another function based on the code of plot.dsa

Best

huipharmacist commented 1 year ago

Is it convenient to send an example