# Simulation of outbreak data

#### Xavier Didelot

#### 2020-06-04

If you want to reproduce exactly the same results as the ones shown in this tutorial, you need to set the seed of your random number generator to zero:

    library(TransPhylo)
    set.seed(0)

A pathogen has an effective within-host population size of \(N_e=100\) and a generation time \(g=1\) day, so that \(N_e g=100/365\) year. The offspring distribution is negative binomial with mean equal to the basic reproduction number \(R=5\). Both the generation time and the sampling time are Gamma distributed with parameters (10,0.1) which has a mean of 1 year. The density of sampling is \(\pi=0.25\). The following commands specify these parameters:

    neg=100/365
    off.r=5
    w.shape=10
    w.scale=0.1
    pi=0.25

We simulate an outbreak that starts in 2005 and which and is observed up to 2008:

    simu <- simulateOutbreak(neg=neg,pi=pi,off.r=off.r,w.shape=w.shape,
                             w.scale=w.scale,dateStartOutbreak=2005,dateT=2008)

This simulation contains both the transmission tree between infected hosts and the within-host phylogenetic tree of each host. This can be visualised as a colored phlogenetic tree, where each host is represented by a unique color:

    plot(simu)

The transmission tree can be extracted and plotted separately from the phylogeny:

    ttree<-extractTTree(simu)
    plot(ttree)

A more detailed plot can be displayed as follows:

    plot(ttree,type='detailed',w.shape,w.scale)

The phylogenetic tree can be extracted and plotted separately from the transmission tree:

    ptree<-extractPTree(simu)
    plot(ptree)

The extracted phylogenetic tree can also be converted into a phylo object from the ape package:

    library(ape)
    p<-phyloFromPTree(ptree)
    plot(p)
    axisPhylo(backward = F)

You can save this tree into a Newick file for further analysis. This is the tree that is used as the starting poit of the tutorial on inference of a transmission tree from a dated phylogeny.

    write.tree(p,'tree.nwk')

The content of this Newick file is:

    write.tree(p,'')

    ## [1] "((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);"

This phylogeny is scaled in years, but time is measured only relatively to the date of the last sample which is at 0 on the x-axis of the figure above. To use this tree again we also need to know exactly when was the last sample taken:

    dateLastSample(simu)

    ## [1] 2007.964