Homec4science
Diffusion hcv_modeling (HCV_b1) Browse

R3744/FP_Scenarios.R
31aa88dcb1eeHCV_b1

Last Change
Raw File

FP_Scenarios.R
View Options

# print("Please enter the interested fibrosis progresssion scenario: 0. Default (age at HCV) 1.FP_Vogl(alcohol and age and gender), 2.FP_Voglmax 3.scenaro_fp_notagedep, 4. scenario_FP_Thein, 5.scenario_maxresFP")
# print("6. scenario_resFP_half, 7. scenario_resFPinOH, 8. scenario_Vogel, 9.scenario_noresidualFP (No residual FP), 10.Default(age at HCV)")
# FP_scenario <- readline(prompt = "FP_scenario:")
FP_scenario <- 1
source("index_block_fs_stage.r") # it identify that in which block and stage we are in
######################################################################################
if (FP_scenario == 0)
{
fp.fun = function(t,lr, bl, history)
{
index.fun = stage_identifier(lr, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
state = index.fun[3]
#######################################
for (state in 31:36)
{
index_fs = 6 # DC
}
######################################
for (state in 37:42)
index_fs = 7 # HCC
####################################
for (state in 43:48)
index_fs = 8 #LT
#####################################################
# identify the state we are in:
ns = 52
m=gems:::auxcounter(ns) # matrix of all transitions (without the terminal states)
if(all(history==0)) {index_trans=1} else {index_trans=max(which(history>0))}
index_state=which(m==index_trans, arr.ind=TRUE)[2] # to get the colum number(arrival state) that corresponds to the last transition
#return in which block is a certain state:
# ref0 = cbind(as.vector(apply(matrix(1:blocks_nb,blocks_nb,1), 1, function(x) rep(x,blocks_size))), 1:(ns-3))
ref1 = cbind(as.vector(apply(matrix(1:blocks_nb,blocks_nb,1), 1, function(x) rep(x,blocks_size))), 1:(statesNum - 22))
ref2 = cbind(c(1,1,1,2,2,2,3,3,3,4,4,4,5,5,5,6,6,6), c(31,37,43,32,38,44,33,39,45,34,40,46,35,41,47,36,42,48))
ref0 = rbind(ref1,ref2)
index_block.fun = function(state, ref=ref0){
ref[which(ref[,2] == state, arr.ind = TRUE),1]
}
index_block=index_block.fun(state = index_state, ref = ref0)
#return in which fibrosis state is a certain state (which level of Fibrosis prog?)
#ref0=cbind(rep(1:blocks_size, blocks_nb), 1:(ns-3))
ref1 = cbind(as.vector(apply(matrix(1:blocks_nb,blocks_nb,1), 1, function(x) rep(x,blocks_size))), 1:(statesNum - 22))
ref2 = cbind(c(1,1,1,2,2,2,3,3,3,4,4,4,5,5,5,6,6,6), c(31,37,43,32,38,44,33,39,45,34,40,46,35,41,47,36,42,48))
ref0 = rbind(ref1,ref2)
index_fs.fun = function(state, ref=ref0){
ref[which(ref[,2] == state, arr.ind = TRUE),1]
}
index_fs=index_fs.fun(state = index_state, ref = ref0)
#######################################################
#####################################
rate_r<-1 # rate ratio
ifelse(index_block == 4 | index_block == 5 | index_block == 6, rate_r <- 0.1, rate_r <- 1) #undetectable HCV (under treatment or SVR/SC) / detectable HCV VL
# at these two blocks we can have the cleaning HCV which can decrease the rate at which fibrosis progressed from F0 to F4
res <- rep(0, length(t))
{
res[bl["Age"]<21] <- { 0.091*rate_r*exp(lr) +0*t }
res[bl["Age"]>=21] <- { 0.105*rate_r*exp(lr) +0*t }
res[bl["Age"]>=31] <- { 0.138*rate_r*exp(lr) +0*t }
res[bl["Age"]>=41] <- { 0.2 *rate_r*exp(lr) +0*t }
res[bl["Age"]>=51] <- { 0.333*rate_r*exp(lr) +0*t }
}
return(res)
}# end of fp.fun
}
##########################################################################################
if (FP_scenario == 1)
{
### FIBROSIS PROGRESSION VOGEL
fp.fun = function(t, bl, history, r)
{
index.fun = stage_identifier(r, bl, history)
index_block = index.fun[1]
index_fs = index.fun[2]
state = index.fun[3]
#######################################
for (state in 31:36)
{
index_fs = 6 # DC
}
######################################
for (state in 37:42)
index_fs = 7 # HCC
####################################
for (state in 43:48)
index_fs = 8 #LT
####################################
#####################################
res<-0 + 0 * t
rr<- 0
{
ifelse(index_block == 1 || index_block == 2 ||index_block == 3 ||index_block == 4 || index_block == 5 , rr<- 1, rr<- 0.1 ) # detectable VL/undetectable VL
{
###################
r = switch (index_fs,
c(1,1.16,1.33),
c(1,1.3,2.22),
c(1,1.3,2.22),
c(1,1.16,4)
)
#r = c(1,1,1.33)
###################
if (index_fs == 1)
{
if (bl["Gender"] == 1) #male
{
res[(t + sum(history) + bl["Age"])< 20] <- { (rr * 0.045 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 20] <- { (rr*0.037 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 30] <-{ (rr*0.027 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 40] <- { (rr * 0.099 * r[bl["Alcohol"]+1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 50] <- { (rr * 0.121 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 60] <- { (rr * 0.138 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 70] <- { (rr * 0.155 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 80] <- { (rr * 0.127 * r[bl["Alcohol"] + 1]) + 0 * t }
}
if (bl["Gender"] == 0) # female
{
res[(t + sum(history) + bl["Age"])< 20] <- { (rr * 0.038 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 20] <- { (rr*0.031 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 30] <-{ (rr*0.022 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 40] <- { (rr * 0.082 * r[bl["Alcohol"]+1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 50] <- { (rr * 0.102 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 60] <- { (rr * 0.115 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 70] <- { (rr * 0.129 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 80] <- { (rr * 0.106 * r[bl["Alcohol"] + 1]) + 0 * t }
}
}#(index_fs == 1)
######################
if (index_fs == 2)
{
if (bl["Gender"] == 1)
{
res[(t + sum(history) + bl["Age"])< 20] <- { (rr * 0.033 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 20] <- { (rr*0.027 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 30] <-{ (rr*0.019 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 40] <- { (rr * 0.072 * r[bl["Alcohol"]+1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 50] <- { (rr * 0.088 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 60] <- { (rr * 0.100 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 70] <- { (rr * 0.112 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 80] <- { (rr * 0.130 * r[bl["Alcohol"] + 1]) + 0 * t }
}
if (bl["Gender"] == 0)
{
res[(t + sum(history) + bl["Age"])< 20] <- { (rr * 0.028 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 20] <- { (rr*0.022 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 30] <-{ (rr*0.016 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 40] <- { (rr * 0.060 * r[bl["Alcohol"]+1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 50] <- { (rr * 0.074 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 60] <- { (rr * 0.083 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 70] <- { (rr * 0.094 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 80] <- { (rr * 0.077 * r[bl["Alcohol"] + 1]) + 0 * t }
}
}#(index_fs == 2)
######################
if (index_fs == 3)
{
if (bl["Gender"] == 1)
{
res[(t + sum(history) + bl["Age"])< 20] <- { (rr * 0.047 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 20] <- { (rr*0.038 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 30] <-{ (rr*0.028 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 40] <- { (rr * 0.102 * r[bl["Alcohol"]+1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 50] <- { (rr * 0.124 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 60] <- { (rr * 0.141 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 70] <- { (rr * 0.159 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 80] <- { (rr * 0.130 * r[bl["Alcohol"] + 1]) + 0 * t }
}
if (bl["Gender"] == 0)
{
res[(t + sum(history) + bl["Age"])< 20] <- { (rr * 0.039 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 20] <- { (rr*0.031 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 30] <-{ (rr*0.023 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 40] <- { (rr * 0.085 * r[bl["Alcohol"]+1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 50] <- { (rr * 0.104 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 60] <- { (rr * 0.118 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 70] <- { (rr * 0.132 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 80] <- { (rr * 0.109 * r[bl["Alcohol"] + 1]) + 0 * t }
}
}#(index_fs == 3)
######################
if (index_fs == 4)
{
if (bl["Gender"] == 1)
{
res[(t + sum(history) + bl["Age"])< 20] <- { (rr * 0.006 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 20] <- { (rr*0.018 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 30] <-{ (rr*0.040 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 40] <- { (rr * 0.063 * r[bl["Alcohol"]+1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 50] <- { (rr * 0.034 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 60] <- { (rr * 0.070 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 70] <- { (rr * 0.136 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 80] <- { (rr * 0.136 * r[bl["Alcohol"] + 1]) + 0 * t }
}
if (bl["Gender"] == 0)
{
res[(t + sum(history) + bl["Age"])< 20] <- { (rr * 0.004 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 20] <- { (rr*0.015 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 30] <-{ (rr*0.033 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 40] <- { (rr * 0.053 * r[bl["Alcohol"]+1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 50] <- { (rr * 0.028 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 60] <- { (rr * 0.059 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 70] <- { (rr * 0.113 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(t + sum(history)+bl["Age"])>= 80] <- { (rr * 0.113 * r[bl["Alcohol"] + 1]) + 0 * t }
}
}#(index_fs == 4)
#######################
}
}
return(res)
} # FIBROSIS PROGRESSION VOGEL
}# end if (FP_scenario == 1)
#######################################################################################
if (FP_scenario == 2)
{
### FIBROSIS PROGRESSION VOGEL
fp.fun = function(t, bl, history, r)
{
# identify the state we are in:
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
index_fs=index_fs.fun(state = index_state, ref = ref0)
res<-0+0*t
if(index_block== 1 || index_block== 2 ||index_block== 3 ||index_block== 4 ||index_block== 5) # detectable VL
{
res[(sum(history))<=9/12]<- (3.6*r)+0*t # r is fixed at 1 in params file
res[(sum(history))>9/12]<- (0.8*r)+0*t
}
else # if undetectable VL
{
res<-0.01*r+0*t
}
return(res)
}
}
#####################################################################################
if (FP_scenario == 3)
{
### FIBROSIS PROGRESSION
fp.fun = function(t, bl, history, r)
{
# identify the state we are in:
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
res<-0+0*t
ifelse (index_block== 1 || index_block== 2 ||index_block== 3 ||index_block== 4 ||index_block== 5 , rr<-1, rr<-0.1) # detectable VL
{
res<- { r[bl["Alcohol"]+1]+0*t } # we fix r=1 (1.33 for OH abusers)
}
return(res)
}
}
#####################################################################################
if (FP_scenario == 4)
{
### FIBROSIS PROGRESSION
fp.fun = function(t, r, bl, history){
# identify the state we are in:
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
# res<-0+0*t
rr<-1
ifelse(index_block== 4 || index_block== 5 ||index_block== 6, rr<-0.1, rr<-1 ) #undetectable HCV (under treatment or SVR/SC) / detectable HCV VL
res <- rep(0, length(t))
{
if (index_fs ==1) res<-rr*0.123+0*t
else if (index_fs ==2) res<-rr*0.113+0*t
else if(index_fs ==3) res<-rr*0.124+0*t
else if (index_fs==4) res<- rr*0.116+0*t
}
return(res)
}
###################################################################################
fp.fun = function(t, bl, history, r){
# identify the state we are in:
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
res<-0+0*t
ifelse(index_block== 1 || index_block== 2 ||index_block== 3 , rr<-1, rr<-0.1) # detectable VL
{
if (index_fs ==1) res<-rr*0.123+0*t
else if (index_fs ==2) res<-rr*0.113+0*t
else if(index_fs ==3) res<-rr*0.124+0*t
else if (index_fs==4) res<- rr*0.116+0*t
}
return(res)
}
}
#####################################################################################
if(FP_scenario == 5)
{
fp.fun = function(t, r, bl, history){
# identify the state we are in:
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
# res<-0+0*t
rr<-1
ifelse(index_block== 3 || index_block== 4, rr<-0.22, rr<-1 ) #undetectable HCV (under treatment or SVR/SC) / detectable HCV VL
res <- rep(0, length(t))
{
res[bl["Age"]<21] <- { 0.091*rr*r +0*t } # we fix r=1 in params
res[bl["Age"]>=21] <- { 0.105*rr*r +0*t }
res[bl["Age"]>=31] <- { 0.138*rr*r +0*t }
res[bl["Age"]>=41] <- { 0.2 *rr*r +0*t }
res[bl["Age"]>=51] <- { 0.333*rr*r +0*t }
}
return(res)
}
}
#####################################################################################
if(FP_scenario == 6)
{
### FIBROSIS PROGRESSION
fp.fun = function(t, bl, history, r)
{
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
res<-0+0*t
rr_fp<-1
ifelse(index_block== 1 || index_block== 2 ||index_block== 3 ||index_block== 5, rr_fp<-1, rr_fp<-0.5) # detectable VL rr_fp =0.22
{
res[bl["Age"]<21] <- { (rr_fp*0.091*r[bl["Alcohol"]+1])+0*t } # we fix r=1 (1.33 for OH abusers)
res[bl["Age"]>=21] <- { (rr_fp*0.105*r[bl["Alcohol"]+1])+0*t }
res[bl["Age"]>=31] <-{ (rr_fp*0.138*r[bl["Alcohol"]+1])+0*t }
res[bl["Age"]>=41] <- { (rr_fp*0.2*r[bl["Alcohol"]+1])+0*t }
res[bl["Age"]>=51] <- { (rr_fp*0.333*r[bl["Alcohol"]+1])+0*t }
}
return(res)
}
}
#################################################################################
if(FP_scenario == 7)
{
### FIBROSIS PROGRESSION
fp.fun = function(t, bl, history, r){
# identify the state we are in:
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
res<-0+0*t
rr_fp<-1
if(index_block== 1 || index_block== 2 ||index_block== 3 ||index_block== 4 ||index_block== 5) rr_fp<-1
else (
ifelse (bl["Alcohol"]==0, rr_fp<-0, rr_fp<-0.1) # detectable VL rr_fp =0.22
)
{
res[bl["Age"]<21] <- { (rr_fp*0.091*r[bl["Alcohol"]+1])+0*t } # we fix r=1 (1.33 for OH abusers)
res[bl["Age"]>=21] <- { (rr_fp*0.105*r[bl["Alcohol"]+1])+0*t }
res[bl["Age"]>=31] <-{ (rr_fp*0.138*r[bl["Alcohol"]+1])+0*t }
res[bl["Age"]>=41] <- { (rr_fp*0.2*r[bl["Alcohol"]+1])+0*t }
res[bl["Age"]>=51] <- { (rr_fp*0.333*r[bl["Alcohol"]+1])+0*t }
}
return(res)
}
}
####################################################################################
if(FP_scenario == 8)
{
### FIBROSIS PROGRESSION MAX as described in VOGEL (2012)
fp.fun = function(t, r, bl, history){
# identify the state we are in:
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
# res<-0+0*t
rr<-1
ifelse(index_block== 3 || index_block== 4 ||index_block== 5 , rr<-0.1, rr<-1 ) #undetectable HCV (under treatment or SVR/SC) / detectable HCV VL
res <- rep(0, length(t))
{
# we fix r=1 in the params file
res<- (rr*0.8*r)+0*t
res[(sum(history))<=0.75]<- (rr*3.6*r)+0*t
}
return(res)
}
}
###################################################################################
if(FP_scenario == 9)
{
#####################################
# No residual FP between F0 and F4
#######################################
fp.fun = function(t, r, bl, history){
# identify the state we are in:
index.fun = stage_identifier(t, r, bl, history)
index_block=index.fun[1]
index_fs = index.fun[2]
# res<-0+0*t
rr<-1
ifelse(index_block== 3 || index_block== 4 ||index_block== 5 , rr<-0, rr<-1 ) #undetectable HCV (under treatment or SVR/SC) / detectable HCV VL
res <- rep(0, length(t))
{
res[bl["Age"]<21] <- { 0.091*rr*r +0*t } # we fix r=1 in params
res[bl["Age"]>=21] <- { 0.105*rr*r +0*t }
res[bl["Age"]>=31] <- { 0.138*rr*r +0*t }
res[bl["Age"]>=41] <- { 0.2 *rr*r +0*t }
res[bl["Age"]>=51] <- { 0.333*rr*r +0*t }
}
return(res)
}
}
#############################################################################
##########################################################################################
if (FP_scenario == 10)
{
### FIBROSIS PROGRESSION VOGEL
fp.fun = function(t, bl, history, r)
{
index.fun = stage_identifier(r, bl, history)
index_block = index.fun[1]
index_fs = index.fun[2]
state = index.fun[3]
#######################################
for (state in 31:36)
{
index_fs = 6 # DC
}
######################################
for (state in 37:42)
index_fs = 7 # HCC
####################################
for (state in 43:48)
index_fs = 8 #LT
#####################################
res<-0 + 0 * t
rr<- 0
{
ifelse(index_block == 1 || index_block == 2 ||index_block == 3 ||index_block == 4 || index_block == 5 , rr<- 1, rr<- 0.1 ) # detectable VL/undetectable VL
{
r = c(1,1.3)
print("chech FP for r")
res[(bl["Age"])< 21] <- { (rr * 0.091 * r[bl["Alcohol"] + 1]) + 0 * t } # we fix r=1 (1.33 for OH abusers)
res[(bl["Age"])>= 21] <- { (rr*0.105 * r[bl["Alcohol"] + 1]) + 0 * t }
res[(bl["Age"])>= 31] <-{ (rr*0.138 * r[bl["Alcohol"] + 1]) + 0*t }
res[(bl["Age"])>= 41] <- { (rr * 0.2 * r[bl["Alcohol"]+1]) + 0 * t }
res[bl["Age"]>= 51] <- { (rr * 0.333 * r[bl["Alcohol"] + 1]) + 0 * t }
# res[(sum(history))<= 2]<- (rr * 3.6 * r) + 0 * t
}
}
return(res)
} # FIBROSIS PROGRESSION VOGEL
}# end if (FP_scenario == 10)
#######################################################################################
c4science · Help