5. ClonalFrameML Integration

The treeWAS R package has also been designed to work with the output of ClonalFrameML. By using the simple read.CFML function, you can convert the output of ClonalFrameML into a form suitable for input into treeWAS.

The Data

For practice using ClonalFrameML, download the ClonalFrameML example data and follow the steps on the wiki to arrive at the example.output dataset. To run the example below, replace the prefix with the path to the example.output dataset on your computer or to another set of output from ClonalFrameML containing:

1. "prefix.labelled_tree.newick"
2. "prefix.ML_sequence.fasta"
3. "prefix.position_cross_reference.txt"

Convert the Data

To convert the data, use the read.CFML function:

prefix <- "./example.output"
dat <- read.CFML(prefix=prefix)

Isolate the elements of the output of read.CFML:

## Required input into treeWAS:
snps <- dat$snps

## Recommended input into treeWAS:
tree <- dat$tree

## Optional input into treeWAS:
n.subs <- dat$n.subs
snps.rec <- dat$snps.rec

Input & Arguments

You can now run treeWAS using this snps matrix, and the tree returned from read.CFML.

Instead of supplying snps.rec to the argument snps.reconstruction, it may be worth (re-)creating the ancestral state reconstruction within treeWAS. This is recommended, provided the dataset you are working with is not very large, in case any inconsistencies exist between the reconstruction performed by ClonalFrameML and that run within treeWAS (which will also be applied to the simulated null dataset).

The dist object can be provided in the n.subs argument, though you may also leave this NULL and have treeWAS recreate n.subs internally, again in case there is any discrepancy between the internal and external methods or reconstruction.

All you would need to run treeWAS with this example dataset is a phenotype. We can make a toy phenotype for this purpose, although no significant findings should result:

set.seed(1)
phen <- sample(c(0,1), nrow(snps), replace=TRUE)
phen <- as.factor(phen)
names(phen) <- rownames(snps)

## Exmine toy phenotype
str(phen)
table(phen)

Running treeWAS

You can now run treeWAS using the converted ClonalFrameML output:

## Example not run... 
out <- treeWAS(snps = snps,
                phen = phen,
                tree = tree,
                n.subs = NULL,
                n.snps.sim = ncol(snps)*10,
                chunk.size = ncol(snps),
                test = c("terminal", "simultaneous", "subsequent"),
                snps.reconstruction = "parsimony",
                snps.sim.reconstruction = "parsimony",
                phen.reconstruction = "parsimony",
                phen.type = NULL,
                na.rm = TRUE,
                p.value = 0.01,
                p.value.correct = "bonf",
                p.value.by = "count",
                dist.dna.model = "JC69",
                plot.tree = FALSE,
                plot.manhattan = TRUE,
                plot.null.dist = TRUE,
                plot.dist = FALSE,
                snps.assoc = NULL, 
                filename.plot = NULL,
                seed = 1)

---