inst/doc/infer.R

In TransPhylo: Inference of Transmission Tree from a Dated Phylogeny

## ----echo=F-------------------------------------------------------------------
knitr::opts_chunk$set(fig.width=6, fig.height=4) 

## -----------------------------------------------------------------------------
library(TransPhylo)
set.seed(0)

## -----------------------------------------------------------------------------
library(ape)
phy<-read.tree(text='((4:1.257652937,(1:1.231048819,5:1.519248672):0.303038892):0.784065883,(3:1.643413444,(6:0.656820028,2:0.007344035611):0.7562780805):1.293120815);')
plot(phy)
axisPhylo(backward = F)

## -----------------------------------------------------------------------------
ptree<-ptreeFromPhylo(phy,dateLastSample=2007.964)
plot(ptree)

## -----------------------------------------------------------------------------
w.shape=10
w.scale=0.1

## -----------------------------------------------------------------------------
dateT=2008

## ----results='hide'-----------------------------------------------------------
res<-inferTTree(ptree,mcmcIterations=1000,w.shape=w.shape,w.scale=w.scale,dateT=dateT)

## -----------------------------------------------------------------------------
plot(res)

## -----------------------------------------------------------------------------
library(coda)
mcmc=convertToCoda(res)
effectiveSize(mcmc)

## -----------------------------------------------------------------------------
med=medTTree(res)
plot(med)

## -----------------------------------------------------------------------------
ttree=extractTTree(med)
plot(ttree,type='detailed',w.shape,w.scale)

## -----------------------------------------------------------------------------
mat=computeMatWIW(res)
lattice::levelplot(mat,xlab='',ylab='')

## -----------------------------------------------------------------------------
mat=computeMatTDist(res)
lattice::levelplot(mat,xlab='',ylab='')

## -----------------------------------------------------------------------------
a=getIncidentCases(res,show.plot = T)

## -----------------------------------------------------------------------------
a=getGenerationTimeDist(res,show.plot = T)

## -----------------------------------------------------------------------------
a=getSamplingTimeDist(res,show.plot = T)

## -----------------------------------------------------------------------------
a=getInfectionTimeDist(res,k=c('1','2'),show.plot = T)

## -----------------------------------------------------------------------------
a=getOffspringDist(res,k=c('1','2'),show.plot = T)