Skip to content

Commit

Permalink
1) added code for model 1 (vaccine model #1)
Browse files Browse the repository at this point in the history
  • Loading branch information
jtm6 committed Jan 9, 2020
1 parent 5ed1161 commit a19ce15
Show file tree
Hide file tree
Showing 6 changed files with 148 additions and 136 deletions.
12 changes: 4 additions & 8 deletions pkg/R/input_parameters_agent_attributes.R
Original file line number Diff line number Diff line change
Expand Up @@ -57,14 +57,10 @@ input_parameters_agent_attributes <-function(){
"virus_part_res_drug", "virus_3_plus_drug_muts",
"virus_3_plus_drug_muts", "Aim3RoundingErrors",
"aim3_mutations_long", "CYP_6_slow",

# -- testing for and treatment of drug resistant viruses (aim 3) --- #
"eligible_2nd_line_ART",
"treated_2nd_line",
"diag_resist_status",
"diag_resist_time",
"last_neg_resist_test",
"time_init_2nd_line",

# -- "new" vaccine code--- #
"phi", "mu",
"sigma",

# -- therapeutic vaccine --- #
"LogSetPoint_genotype", "vacc_status_at_inf",
Expand Down
4 changes: 3 additions & 1 deletion pkg/R/transmission_hill_fxn.R
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,9 @@ result <- result * susceptibility # New for AgeAndSPVL work
result[which(acute_phase_status==1)] <- result[which(acute_phase_status==1)] * dat$param$trans_RR_acute_phase
result[which(sti_status_sus==1)] <- result[which(sti_status_sus==1)] * dat$param$trans_RR_STI
result[which(condom_use==1)] <- result[which(condom_use==1)] * dat$param$trans_RR_uses_condoms
result[which(vacc_sens_sus==1)] <- result[which(vacc_sens_sus==1)] *dat$param$trans_RR_vaccine

#for vaccination branch, vaccination dynamics occurs elsewhere, comment out vaccine effects
#result[which(vacc_sens_sus==1)] <- result[which(vacc_sens_sus==1)] *dat$param$trans_RR_vaccine

# Age effect
temp_age_xb <- dat$param$trans_RR_age^((age_vec_sus-dat$param$trans_base_age)/10)
Expand Down
3 changes: 2 additions & 1 deletion pkg/R/transmission_hughes_and_exp.R
Original file line number Diff line number Diff line change
Expand Up @@ -57,8 +57,9 @@ transmission_hughes_and_exp<-function(dat, acute_phase_status, sti_status_sus, c
# Age effect
temp_age_xb <- (age_vec_sus - dat$param$trans_base_age) * log(dat$param$trans_RR_age)/10
xB <- xB + temp_age_xb
#for vaccination branch, vaccination dynamics occurs elsewhere, comment out vaccine effects
#vaccine effect
xB <- xB + vacc_sens_sus * log(dat$param$trans_RR_vaccine)
#xB <- xB + vacc_sens_sus * log(dat$param$trans_RR_vaccine)
exp_xB = exp(xB)

result <- 1 - ((1-dat$param$trans_lambda)^(1*exp_xB))
Expand Down
49 changes: 13 additions & 36 deletions pkg/R/transmission_main_module.R
Original file line number Diff line number Diff line change
Expand Up @@ -86,24 +86,6 @@ transmission_main_module <- function(dat,at)
if(length(acute_phase_index)>0){ acute_phase_status[acute_phase_index] <- 1 }
#sti status of susceptible (relevant for both msm and hetero)
sti_status_sus <- dat$pop$sti_status[sus_id]
#on prep

#--------- preventative vaccine dynamics -----------------
#note: multi-efficacy-dynanmics are after calculation of transmission probs ~ line

#set vaccine effect to zero for all agents, but then fill in if necessary
vacc_sens_sus <- rep(0,nrow(dat$discord_coital_df))

#if sus is vaccinated and infected with sensitive virus
if(at >= dat$param$start_vacc_campaign[1] & dat$param$preventative_campaign==T){
vacc_ix <- which(dat$pop$vaccinated[sus_id]==1 &
dat$pop$virus_sens_vacc[inf_id]==1)
if(length(vacc_ix)>0){vacc_sens_sus[vacc_ix] <- 1}
}


#------end of preventative vaccine--------------------



if(dat$param$model_sex=="msm"){
Expand Down Expand Up @@ -143,36 +125,31 @@ transmission_main_module <- function(dat,at)
} #end of hetero variables

if(dat$param$transmission_model=="hill"){
trans_probs <- transmission_hill_fxn(dat, acute_phase_status, sti_status_sus, condom_use,
trans_probs_raw <- transmission_hill_fxn(dat, acute_phase_status, sti_status_sus, condom_use,
msm_sus_receptive_status, msm_sus_insert_status,
msm_circum_status_insert_sus, age_vec_sus,
hetero_sus_female_vi, hetero_sus_female_ai,
hetero_sus_male_ai, hetero_sus_male_circum_status,
V_inf, vacc_sens_sus, susceptibility)
V_inf, vacc_sens_sus=rep(0,nrow(dat$discord_coital_df)),
susceptibility)
}else{
trans_probs <- transmission_hughes_and_exp(dat, acute_phase_status, sti_status_sus, condom_use,
trans_probs_raw <- transmission_hughes_and_exp(dat, acute_phase_status, sti_status_sus, condom_use,
msm_sus_receptive_status, msm_sus_insert_status,
msm_circum_status_insert_sus, age_vec_sus,
hetero_sus_female_vi, hetero_sus_female_ai,
hetero_sus_male_ai, hetero_sus_male_circum_status,
logV_inf, vacc_sens_sus, susceptibility)
logV_inf, vacc_sens_sus=rep(0,nrow(dat$discord_coital_df)),
susceptibility)
}
###################################
# multi-efficacy dynamics
#multi efficacy model
#find which suscept. agent is vaccinated, get their vaccine efficacy
# and fill in "vacc_efficacies" vector (all other values of this vector are zero)
if(at >= dat$param$start_vacc_campaign[1] & dat$param$vacc_multi_eff==T){
vacc_efficacies <- rep(0,nrow(dat$discord_coital_df))
vacc_ix <- which(dat$pop$vaccinated[sus_id]==1)
if(length(vacc_ix)>0){
#for vaccinated and susceptible agents, get partner's VE
inf_ix= inf_id[vacc_ix]
vacc_efficacies[vacc_ix] <- dat$pop$vacc_eff[inf_ix]
trans_probs <- trans_probs*(1-vacc_efficacies[vacc_ix])
}
}
#Vaccine dynamics

#calculate theta based on user-specified vaccine model (1,2,3,...)
temp_fxn<- paste("caclulate_theta",dat$param$vaccine_model_id,sep="")
theta <- do.call(temp_fxn,list(dat))

#adjust raw transmission probabilities
trans_probs <- trans_probs_raw*(1-theta)

###################################
# fill in discord_coital_df
Expand Down
90 changes: 0 additions & 90 deletions pkg/R/vaccination.R

This file was deleted.

126 changes: 126 additions & 0 deletions pkg/R/vaccination_model_1.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,126 @@
#file contains: phi1
# update_mu1
# covariates1
# marks1
# draw_m1
# caclulate_theta1

#################

phi1 <- function(dat, at) {
#1) Once vaccination campaign starts, fxn changes vaccination status (attribute "phi") from
# NA to "1" to eligibe/selected agents.
#2) After vacc. campaign has started, will change status of vaccinated agent's whose
#efficacy has waned to "0" (stochastic draw)


#skip vaccination routine if vacc. campaign hasn't started
if(at < dat$param$start_vacc_campaign[1]) {return(dat)}

# off/on for already vaccinated
if(at > dat$param$start_vacc_campaign[1] ) {
vacc_ix <- which(dat$pop$phi == 1)
dat$pop$phi[vacc_ix] <- rbinom(length(vacc_ix), 1, 1 - (1/dat$param$vacc_eff_duration))
}


#if vaccination only occurs at discrete intervals, skip rest of routine if vaccination
#campaign has started but does not occur in this time step
if(!is.element(at,dat$param$start_vacc_campaign)) {return(dat)}
#------------------------------------------------

# If vaccine is targeted to attribute groups

#if designated vacc. level reached (percent of pop vaccianted), don't vacc anymore
if(length(which(dat$pop$phi == 1 ))/length(which(dat$pop$Status>=0)) > dat$param$max_perc_vaccinated){return(dat)}


# Eligible_patients: eligible for care, not vaccinated, not infected
eligible_index <- which(dat$pop$Status == 0 &
(dat$pop$phi == 0 | is.na(dat$pop$phi)) &
dat$pop$eligible_care == 1)

if(length(eligible_index) == 0) {return(dat)} #if no agents are eligible

no_vaccinated <- sum(rbinom(length(which(dat$pop$Status>=0)), 1, dat$param$perc_vaccinated)) #denominator is total population alive
if(no_vaccinated == 0) {return(dat)}


if(no_vaccinated <length(eligible_index)){
vaccinated_index <- sample(eligible_index, no_vaccinated)
}else{
vaccinated_index <- eligible_index
#if the %coverage in total population alive exceeds #eligible, vaccinate all eligible
}

dat$pop$phi[vaccinated_index] <- 1
dat$pop$vacc_init_time[vaccinated_index] <- at

return(dat)
}
#######################

update_mu1 <- function (dat,at){
#code to initialize/update mu,sigma for each agent or initialize for new agents to model
# for model 1, this fxn just returns dat object without modification (no mu,sigma to update)

#note: "m" and "sigma" are agent attributes (see file/fxn "input_parameters_agent_attributes)
#which means that the default value of the attribute is "NA" until modified

return(dat)
}
#######################

covariates1 <- function(dat,at){
#initializupdate baseline covariates necessary for vaccine model's "calculate_theta" fxn
#that are not already part of baseline evonete/
#for model 1, no covariates so just return "dat$covariates" object unmodificed
if(at==2){
dat$covariates <- NA
}

return(dat)
}
######################

marks1 <- function(dat,at){
#at start of model, initialize "marks" for each agent
#after start, for given agent specific mu/sigma in dat$mu and dat$sigma,
#calculate updated marks,stored in dat$m

#for model1, marks for each agent is static, that is, agent is either infected or not
#thus, for model1, marks is just the "Status" attribute (0/1)

dat$pop$m <- dat$pop$Status
return(dat)

}
######################
draw_m1 <- function(dat,at){
# calculate marks, based on most current values of mu and sigma

#for model 1, marks don't change, it is whether infector's virus is sensitive
#to infectee's vaccine (0/1)
m <- dat$pop$virus_sens_vacc[dat$infector_id]
return(m)
}
#######################

caclulate_theta1 <- function(dat){
# draw updated marks (of virus)
# then caclulate theta based on dat$m, dat$phi, dat$covariates
# output theta to adjust raw trans probs

#for model 1, mark is whether virus of infector is sensitive to vaccine of
# infectee

m <- draw_m1(dat,at)
#theta is vaccination status (0/1) *
# virus susceptibility to vaccine (0/1) *
# percent reduction due to vaccine in trans. prob.

theta <- dat$pop$phi[dat$susceptible_id]*m*dat$param$vacc_trans_prob_decrease
return(theta)
}


0 comments on commit a19ce15

Please sign in to comment.