#######################################################
#######################################################
####### #######
####### CNA and SNA input #######
####### #######
#######################################################
#######################################################
Library (Canopy)
Data ("MDA231")
ProjectName = mda231$projectname # # Name of Project
R = Mda231$r; R # # Mutant allele read depth (for SNAs)
X = mda231$x; X # # Total Depth (for SNAs)
WM = MDA231$WM; WM # # observed major copy number (for CNA regions)
Wm = MDA231$WM; Wm # # observed minor copy number (for CNA regions)
Epsilonm = mda231$epsilonm # # Standard deviation of the WM, pre-fixed here
Epsilonm = mda231$epsilonm # # Standard deviation of the Wm, pre-fixed here
# # Whether CNA Regions Harbor specific CNAs (only needed for overlapping cnas)
C = mda231$c; C
Y = mda231$y; Y # # Whether SNAs is affected by CNAs
#######################################################
#######################################################
####### #######
####### MCMC Sampling #######
####### #######
#######################################################
#######################################################
K = 3:6 # of Subclones
Numchain = # Number of chains with random initiations
Sampchain = Canopy.sample (R = r, x = x, wm = WM, WM = WM, Epsilonm = Epsilonm,
Epsilonm = epsilonm, c = c, y = y, k = k,
Numchain = Numchain, Simrun = 100000, Writeskip = 200,
ProjectName = ProjectName, Cell.line = TRUE,
Plot.likelihood = TRUE)
Save.image (file = paste (ProjectName, ' _postmcmc_image.rda ', sep= '),
Compress = ' XZ ')
#######################################################
#######################################################
####### #######
####### BIC to determine number of Subclones #######
####### #######
#######################################################
#######################################################
Library (Canopy)
Projectname= ' MDA231 '
Load (paste (ProjectName, ' _postmcmc_image.rda ', sep= '))
Burnin = 100
Thin = 10
# If PDF = TRUE, a PDF would be generated.
BIC = canopy. BIC (Sampchain = sampchain, ProjectName = projectname, k = k,
Numchain = numchain, Burnin = burnin, thin = thin, pdf = TRUE)
OPTK = K[which.max (BIC)]
#######################################################
#######################################################
####### #######
####### Posterior tree Evaluation #######
####### #######
#######################################################
#######################################################
Post = canopy.post (Sampchain = sampchain, ProjectName = projectname, k = k,
Numchain = numchain, Burnin = burnin, thin = thin,
OPTK = OPTK, c = c, Post.config.cutoff = 0.05)
Samptreethin = post[[1]] # List of all post-burnin and thinning trees
Samptreethin.lik = post[[2] "# Likelihoods of Trees in Samptree
Config = post[[3]]
Config.summary = Post[[4]]
Print (config.summary)
# First Column:tree configuration
# Second Column:posterior configuration probability in the entire tree space
# Third Column:posterior configuration likelihood in the subtree space
# note:if modes of posterior probabilities aren ' t obvious, run sampling longer.
#######################################################
#######################################################
####### #######
####### Tree output and plot #######
####### #######
#######################################################
#######################################################
# Choose the configuration with the highest posterior likelihood
CONFIG.I = Config.summary[which.max (config.summary[,3]), 1]
Cat (' Configuration ', CONFIG.I, ' has the highest posterior likelihood.\n ')
Output.tree = Canopy.output (post, CONFIG.I, C)
Pdf.name = Paste (ProjectName, ' _config_highest_likelihood.pdf ', sep= ')
Canopy.plottree (output.tree, pdf = TRUE, pdf.name = pdf.name)
# Plot posterior tree with second configuration
Output.tree = Canopy.output (post, 1, C)
Canopy.plottree (Output.tree, pdf=true, pdf.name = paste (ProjectName, ' _second_config.pdf ', Sep = '))
A history of "evolutionary tree" of tumors based on exon data from cancerous cells