## ----setup, include = FALSE--------------------------------------------------- knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) ## ----getHostData, eval = FALSE------------------------------------------------ # getHostData(nosoi.output, what, pop) ## ----getTableHosts, eval = FALSE---------------------------------------------- # getTableHosts(nosoi.output, pop) ## ----getTableState, eval = FALSE---------------------------------------------- # getTableState(nosoi.output, pop) ## ----nosoiSummary, eval = FALSE----------------------------------------------- # summary(nosoi.output) ## ----simDynamics, message=FALSE----------------------------------------------- library(nosoi) t_incub_fct <- function(x){rnorm(x,mean = 5,sd=1)} p_max_fct <- function(x){rbeta(x,shape1 = 5,shape2=2)} p_Move_fct <- function(t){return(0.1)} p_Exit_fct <- function(t){return(0.05)} proba <- function(t,p_max,t_incub){ if(t <= t_incub){p=0} if(t >= t_incub){p=p_max} return(p) } time_contact <- function(t, current.in, host.count){ temp.val = 30 - host.count if(temp.val <= 0) { return(0) } if(temp.val >= 0) { if(current.in=="A"){ return(round((temp.val/30)*rnorm(1, 3, 1), 0))} if(current.in=="B"){return(0)} if(current.in=="C"){ return(round((temp.val/30)*rnorm(1, 6, 1), 0))} } } transition.matrix = matrix(c(0,0.2,0.4,0.5,0,0.6,0.5,0.8,0),nrow = 3, ncol = 3,dimnames=list(c("A","B","C"),c("A","B","C"))) set.seed(1050) test.nosoiA <- nosoiSim(type="single", popStructure="discrete", length=100, max.infected=200, init.individuals=1, init.structure="A", structure.matrix=transition.matrix, pMove=p_Move_fct, param.pMove=NA, diff.nContact=TRUE, hostCount.nContact=TRUE, nContact=time_contact, param.nContact=NA, pTrans = proba, param.pTrans = list(p_max=p_max_fct, t_incub=t_incub_fct), pExit=p_Exit_fct, param.pExit=NA ) ## ----figureDynamics1, eval=FALSE, message=FALSE------------------------------- # library(ggplot2) # cumulative.table <- getCumulative(test.nosoiA) # dynamics.table <- getDynamic(test.nosoiA) # # ggplot(data=cumulative.table, aes(x=t, y=Count)) + geom_line() + theme_minimal() + labs(x="Time (t)",y="Cumulative count of infected hosts") ## ----figureDynamics1_actual, echo=FALSE, message=FALSE------------------------ if (!requireNamespace("ggplot2", quietly = TRUE)) { message("Package 'ggplot2' is needed for plotting this figure.") } else { library(ggplot2) cumulative.table <- getCumulative(test.nosoiA) dynamics.table <- getDynamic(test.nosoiA) ggplot(data=cumulative.table, aes(x=t, y=Count)) + geom_line() + theme_minimal() + labs(x="Time (t)",y="Cumulative count of infected hosts") } ## ----figureDynamics2, eval=FALSE, message=FALSE------------------------------- # ggplot(data=dynamics.table, aes(x=t, y=Count, color=state)) + geom_line() + theme_minimal() + labs(x="Time (t)",y="Number of active infected hosts") ## ----figureDynamics2_actual, echo=FALSE, message=FALSE------------------------ if (!requireNamespace("ggplot2", quietly = TRUE)) { message("Package 'ggplot2' is needed for plotting this figure.") } else { ggplot(data=dynamics.table, aes(x=t, y=Count, color=state)) + geom_line() + theme_minimal() + labs(x="Time (t)",y="Number of active infected hosts") } ## ----simDynamics-bis, eval=FALSE, message=FALSE------------------------------- # library(nosoi) # t_incub_fct <- function(x){rnorm(x,mean = 5,sd=1)} # p_max_fct <- function(x){rbeta(x,shape1 = 5,shape2=2)} # p_Move_fct <- function(t){return(0.1)} # # p_Exit_fct <- function(t){return(0.05)} # # proba <- function(t,p_max,t_incub){ # if(t <= t_incub){p=0} # if(t >= t_incub){p=p_max} # return(p) # } # # time_contact <- function(t, current.in, host.count){ # # temp.val = 30 - host.count # # if(temp.val <= 0) { # return(0) # } # if(temp.val >= 0) { # if(current.in=="A"){ # return(round((temp.val/30)*rnorm(1, 3, 1), 0))} # if(current.in=="B"){return(0)} # if(current.in=="C"){ # return(round((temp.val/30)*rnorm(1, 6, 1), 0))} # } # } # # transition.matrix = matrix(c(0,0.2,0.4,0.5,0,0.6,0.5,0.8,0),nrow = 3, ncol = 3,dimnames=list(c("A","B","C"),c("A","B","C"))) # # set.seed(1050) # test.nosoiA <- nosoiSim(type="single", popStructure="discrete", # length=100, # max.infected=200, # init.individuals=1, # init.structure="A", # structure.matrix=transition.matrix, # pMove=p_Move_fct, # param.pMove=NA, # diff.nContact=TRUE, # hostCount.nContact=TRUE, # nContact=time_contact, # param.nContact=NA, # pTrans = proba, # param.pTrans = list(p_max=p_max_fct, # t_incub=t_incub_fct), # pExit=p_Exit_fct, # param.pExit=NA # ) ## ----R0-1, message=FALSE------------------------------------------------------ getR0(test.nosoiA) ## ----figureR0, eval=FALSE, message=FALSE-------------------------------------- # data = data.frame(R0=getR0(test.nosoiA)$R0.dist) # ggplot(data=data, aes(x=R0)) + geom_histogram() + theme_minimal() ## ----figureR0_actual, echo=FALSE, message=FALSE------------------------------- data = data.frame(R0=getR0(test.nosoiA)$R0.dist) if (!requireNamespace("ggplot2", quietly = TRUE)) { message("Package 'ggplot2' is needed for plotting this figure.") } else { ggplot(data=data, aes(x=R0)) + geom_histogram() + theme_minimal() } ## ----simDynamics-ter, eval=FALSE, message=FALSE------------------------------- # library(nosoi) # t_incub_fct <- function(x){rnorm(x,mean = 5,sd=1)} # p_max_fct <- function(x){rbeta(x,shape1 = 5,shape2=2)} # p_Move_fct <- function(t){return(0.1)} # # p_Exit_fct <- function(t){return(0.05)} # # proba <- function(t,p_max,t_incub){ # if(t <= t_incub){p=0} # if(t >= t_incub){p=p_max} # return(p) # } # # time_contact <- function(t, current.in, host.count){ # # temp.val = 30 - host.count # # if(temp.val <= 0) { # return(0) # } # if(temp.val >= 0) { # if(current.in=="A"){ # return(round((temp.val/30)*rnorm(1, 3, 1), 0))} # if(current.in=="B"){return(0)} # if(current.in=="C"){ # return(round((temp.val/30)*rnorm(1, 6, 1), 0))} # } # } # # transition.matrix <- matrix(c(0,0.2,0.4,0.5,0,0.6,0.5,0.8,0),nrow = 3, ncol = 3,dimnames=list(c("A","B","C"),c("A","B","C"))) # # set.seed(1050) # test.nosoiA <- nosoiSim(type="single", popStructure="discrete", # length=100, # max.infected=200, # init.individuals=1, # init.structure="A", # structure.matrix=transition.matrix, # pMove=p_Move_fct, # param.pMove=NA, # diff.nContact=TRUE, # hostCount.nContact=TRUE, # nContact=time_contact, # param.nContact=NA, # pTrans = proba, # param.pTrans = list(p_max=p_max_fct, # t_incub=t_incub_fct), # pExit=p_Exit_fct, # param.pExit=NA # ) ## ----treeA, eval=FALSE, message=FALSE----------------------------------------- # test.nosoiA.tree <- getTransmissionTree(test.nosoiA) # # if (requireNamespace("ggtree", quietly = TRUE) && # requireNamespace("ggplot2", quietly = TRUE)) { # library(ggtree) # library(ggplot2) # ggtree(test.nosoiA.tree, color = "gray30") + # geom_nodepoint(aes(color=state)) + # geom_tippoint(aes(color=state)) + # theme_tree2() + xlab("Time (t)") + # theme(legend.position = c(0,0.8), # legend.title = element_blank(), # legend.key = element_blank()) # } else { # message("Packages 'ggtree' and 'ggplot2' are needed for plotting") # } ## ----treeA_actual, echo=FALSE, message=FALSE---------------------------------- if (!(requireNamespace("ape", quietly = TRUE) && requireNamespace("tidytree", quietly = TRUE) && requireNamespace("treeio", quietly = TRUE))) { message("Packages 'ape', 'tidytree' and 'treeio' are needed for transmission tree generation.") } else { test.nosoiA.tree <- getTransmissionTree(test.nosoiA) } if (!(requireNamespace("ggplot2", quietly = TRUE) && requireNamespace("ggtree", quietly = TRUE) && (utils::packageVersion("ggtree") < "3.1.3"))) { message("Packages 'ggplot2' and 'ggtree (> 3.1.3)' are needed for plotting this figure.") } else { library(ggplot2) library(ggtree) ggtree(test.nosoiA.tree, color = "gray30") + geom_nodepoint(aes(color=state)) + geom_tippoint(aes(color=state)) + theme_tree2() + xlab("Time (t)") + theme(legend.position = c(0,0.8), legend.title = element_blank(), legend.key = element_blank()) } ## ----tree-sample1, echo=FALSE, message=FALSE---------------------------------- if (!(requireNamespace("ape", quietly = TRUE) && requireNamespace("tidytree", quietly = TRUE) && requireNamespace("treeio", quietly = TRUE))) { message("Packages 'ape', 'tidytree' and 'treeio' are needed for transmission tree generation.") } else { library(dplyr) table.hosts.testA <- getTableHosts(test.nosoiA) setkey(table.hosts.testA, "hosts.ID") set.seed(30082019) sampled.hosts <- sample(table.hosts.testA$hosts.ID, 20) samples.data.table <- table.hosts.testA[sampled.hosts] %>% group_by(hosts.ID) %>% mutate(Sampled.time = round(runif(1, min=inf.time,max=out.time), digits = 0), labels = paste(hosts.ID,"sampled",sep="-")) samples.data.table <- samples.data.table[c("hosts.ID","Sampled.time","labels")] colnames(samples.data.table) <- c("hosts", "times", "labels") samples.data.table } ## ----tree-sample2, eval=FALSE, message=FALSE---------------------------------- # test.nosoiA.tree.sampled <- sampleTransmissionTree(test.nosoiA, test.nosoiA.tree, samples.data.table) ## ----tree-sample2_actual, echo=FALSE, message=FALSE--------------------------- if (!(requireNamespace("ape", quietly = TRUE) && requireNamespace("tidytree", quietly = TRUE) && requireNamespace("treeio", quietly = TRUE))) { message("Packages 'ape', 'tidytree' and 'treeio' are needed for transmission tree generation.") } else { test.nosoiA.tree.sampled <- sampleTransmissionTree(test.nosoiA, test.nosoiA.tree, samples.data.table) } ## ----treeB, eval = FALSE, message=FALSE--------------------------------------- # if (!requireNamespace("ggtree", quietly = TRUE)) { # ggtree(test.nosoiA.tree.sampled, color = "gray30") + geom_nodepoint(aes(color=state)) + geom_tippoint(aes(color=state)) + geom_tiplab(aes(label=host)) + # theme_tree2() + xlab("Time (t)") + theme(legend.position = c(0,0.8), # legend.title = element_blank(), # legend.key = element_blank()) # } ## ----treeB_actual, echo=FALSE, message=FALSE---------------------------------- if (!(requireNamespace("ggtree", quietly = TRUE) && (utils::packageVersion("ggtree") < "3.1.3"))) { message("Package 'ggtree' is needed for plotting this figure.") } else { ggtree(test.nosoiA.tree.sampled, color = "gray30") + geom_nodepoint(aes(color=state)) + geom_tippoint(aes(color=state)) + geom_tiplab(aes(label=host)) + theme_tree2() + xlab("Time (t)") + theme(legend.position = c(0,0.8), legend.title = element_blank(), legend.key = element_blank()) } ## ----tree-sample3, echo=FALSE, message=FALSE---------------------------------- if (!(requireNamespace("ape", quietly = TRUE) && requireNamespace("tidytree", quietly = TRUE) && requireNamespace("treeio", quietly = TRUE))) { message("Packages 'ape', 'tidytree' and 'treeio' are needed for transmission tree generation.") } else { set.seed(5905950) sampled.hosts2 = subset(table.hosts.testA, active == 0)$hosts.ID sampled.hosts <- sample(sampled.hosts2, 20) sampled.hosts } ## ----tree-sample4, eval=FALSE, message=FALSE---------------------------------- # test.nosoiA.tree.sampled.exiting <- sampleTransmissionTreeFromExiting(test.nosoiA.tree, sampled.hosts) ## ----tree-sample4_actual, echo=FALSE, message=FALSE--------------------------- if (!(requireNamespace("ape", quietly = TRUE) && requireNamespace("tidytree", quietly = TRUE) && requireNamespace("treeio", quietly = TRUE))) { message("Packages 'ape', 'tidytree' and 'treeio' are needed for transmission tree generation.") } else { test.nosoiA.tree.sampled.exiting <- sampleTransmissionTreeFromExiting(test.nosoiA.tree, sampled.hosts) } ## ----treeC, eval=FALSE, message=FALSE----------------------------------------- # if (requireNamespace("ggtree", quietly = TRUE)) { # ggtree(test.nosoiA.tree.sampled.exiting, color = "gray30") + geom_nodepoint(aes(color=state)) + geom_tippoint(aes(color=state)) + geom_tiplab(aes(label=host)) + # theme_tree2() + xlab("Time (t)") + theme(legend.position = c(0,0.8), # legend.title = element_blank(), # legend.key = element_blank()) # } else { # message("Package 'ggtree' required for plotting") # } ## ----treeC_actual, echo=FALSE, message=FALSE---------------------------------- if (!(requireNamespace("ape", quietly = TRUE) && requireNamespace("tidytree", quietly = TRUE) && requireNamespace("treeio", quietly = TRUE) && requireNamespace("ggtree", quietly = TRUE) && (utils::packageVersion("ggtree") < "3.1.3"))) { message("Packages 'ape', 'tidytree', 'treeio' and 'ggtree' are needed for transmission tree generation.") } else { ggtree(test.nosoiA.tree.sampled.exiting, color = "gray30") + geom_nodepoint(aes(color=state)) + geom_tippoint(aes(color=state)) + geom_tiplab(aes(label=host)) + theme_tree2() + xlab("Time (t)") + theme(legend.position = c(0,0.8), legend.title = element_blank(), legend.key = element_blank()) } ## ----writeTree, eval=FALSE, message=FALSE------------------------------------- # treeio::write.beast(test.nosoiA.tree.sampled.exiting)