Main
Inferring the clonal structures by the Clonevol (https://github.com/hdng/clonevol).
Inferring the clonal structures by the Clonevol (https://github.com/hdng/clonevol).
library(clonevol)
library(ggrepel)
## Loading required package: ggplot2
library(tidyverse)
## Warning in system("timedatectl", intern = TRUE): running command 'timedatectl'
## had status 1
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✓ tibble 3.1.5 ✓ dplyr 1.0.7
## ✓ tidyr 1.1.4 ✓ stringr 1.4.0
## ✓ readr 2.0.2 ✓ forcats 0.5.1
## ✓ purrr 0.3.4
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
mutdata.vaf1 = readRDS("data/mutdata.vaf1.rds")
new.label = readRDS(file = "data/new.label.rds")
mutdata = mutdata.vaf1 %>%
left_join(new.label)
## Joining, by = c("mutid", "sciClone", "kmeans")
mutdata[, 2:34] = mutdata[, 2:34]*100
mutdata$cluster = mutdata$Clone2
#Visualizing the variant clusters
set.colors = c("#C6C6C6", "#6FDCBF","#5AA36E","#E99692","#B4D985","#EA7D23","#E53CFB","#4B85B7","#8439FA","#BD8736","#B3B371","#A7C7DE","#EE97FC","#57C222","#BFABD0","#44589B", "#794C18", RColorBrewer::brewer.pal(n = 10, name = "Paired") )
clone.colors = set.colors[ 1:length(unique(mutdata$cluster))]
cluster.col.name = "cluster"
head(mutdata)
Labeling the driver mutations.
#Stomach drivers
GC01 = readRDS("data/GC01.rds")
GC01 = GC01 %>% mutate(Chromosome = str_c("chr", Chromosome) ) %>%
mutate(mutid = str_c(Chromosome, Start_Position, Reference_Allele, sep = ":"))
drivers = read.delim("data/IntOGen-Drivers-Cancer_Genes.tsv", header = T) %>%
filter(CANCER_TYPE %in% c("BRCA","COREAD","LUAD", "LUSC") ) %>%
pull(SYMBOL) %>% unique()
mutinfo = GC01 %>%
select(Chromosome, Start_Position, End_Position, Hugo_Symbol, Variant_Classification, Protein_Change, Reference_Allele, Tumor_Seq_Allele2 ) %>%
unique.data.frame() %>%
mutate(#Chromosome = str_c("chr", Chromosome),
Protein_Change = ifelse(Variant_Classification == "Splice_Site", "splice", Protein_Change),
Protein_Change = ifelse(Variant_Classification == "3'UTR", "3UTR", Protein_Change),
Protein_Change = ifelse(Variant_Classification == "5'UTR", "5UTR", Protein_Change),
gene_site = str_c(Hugo_Symbol, "_" ,Protein_Change),
is.driver = ifelse( (Hugo_Symbol %in% drivers) & ( !Variant_Classification %in% c("Silent", "3'Flank", "IGR", "Intron", "RNA") ) , TRUE, FALSE) ) %>%
mutate( mutid = str_c(Chromosome , Start_Position , Reference_Allele , sep = ":") )
mutdata = right_join(mutinfo, mutdata ) %>%
arrange(cluster) %>%
as.data.frame()
## Joining, by = "mutid"
#add ref and var columns for merge data.
sampleNames = c( paste0("Breast_", 1:5), paste0("Lung_", 1:5),
paste0("Coad_", 1:5), paste0("OveryLM_", 1:5),
paste0("OveryRM_", 1:6), paste0("UterusM_", c(1:7)) )
var.data = round(mutdata[, sampleNames])
ref.data = 100 - var.data
colnames(var.data) = paste0(colnames(var.data), ".var")
colnames(ref.data) = paste0(colnames(ref.data), ".ref")
mutdata = cbind(mutdata, var.data, ref.data)
message("<2> Analysis using clonevol")
## <2> Analysis using clonevol
output = "Met2" #Output folder names
#saveRDS(mutinfo, file = sprintf("%s/%s.mutinfo.rds", output))
sampleNames = c( paste0("Coad_", 1:5), paste0("OveryLM_", 1:5),
paste0("OveryRM_", 1:6),paste0("UterusM_", c(1, 3)) )
vaf.col.names = sampleNames
show.ccf = F
pp <- plot.variant.clusters(mutdata,
show.cluster.size = F,
show.cluster.label= F,
cluster.col.name = cluster.col.name,
vaf.col.names = vaf.col.names,
violin = FALSE,
box = FALSE,
jitter = TRUE,
jitter.shape = 1,
variant.class.col.name =cluster.col.name,
#vaf.limits =
jitter.color = clone.colors,
jitter.size = 1.2,
jitter.alpha = 1,
jitter.width = 0.2,
jitter.center.method = "median",
jitter.center.size = 1,
jitter.center.color = "darkgray",
jitter.center.display.value = "none",
display.plot = F,
horizontal = T,
order.by.total.vaf = F,
highlight = 'is.driver',
highlight.shape = 21,
highlight.color = 'blue',
highlight.fill.color = 'green',
highlight.size = 2.5,
highlight.note.col.name = NULL,
highlight.note.size = 2,
highlight.note.color = "blue",
highlight.note.angle = 0,
founding.cluster = 1,
ccf = show.ccf
)
##
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
##
## combine
#if(show.ccf){max.limits = 180 }else{max.limits = 80}
max.limits = 110
pp[[1]] = pp[[1]] +
theme(axis.title.y = element_blank()) +
scale_y_continuous(position = "right", limits = c(0, max.limits)) +
theme(axis.text.x = element_text(angle = 270, size = 8) )
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
for(ii in 2:length(pp)){
pp[[ii ]] = pp[[ii]] +
scale_y_continuous(position = "right", limits = c(0, max.limits)) +
theme(axis.text.x = element_text(angle = 270, size = 8 ),
axis.text.y = element_blank()
)
}
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
#annotate the driver genes.
if(any(mutdata$is.driver) ){
labels = pp[[1]]$data
#select colors for annotation.
clone.colors.sel =unique( clone.colors[ labels$cluster[labels$is.driver]])
pp1 = labels %>% filter(is.driver) %>%
mutate(cluster_1 = factor(cluster)) %>%
ggplot(aes(y = cluster, x = 1, label = gene_site)) + theme_classic() +
geom_point(aes( color = cluster_1), size = 3)+
geom_text_repel(
nudge_x = 0.15,
direction = "y",
hjust = 0,
segment.size = 0.2,
size = 3,
) +
ylim(0, length(clone.colors) +1 ) +
xlim(1, 0.8) +
scale_color_manual( values= clone.colors.sel, guide="none") +
theme(
axis.line = element_blank(),
axis.ticks = element_blank(),
axis.text = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
plot.title = element_text(hjust = 0.5)
)
pp[[ length(pp) +1 ]] = pp1
}
pdf( sprintf( "%s/%s.box.pdf", output, output), width = 2*length(pp), height = 4)
ggpubr::ggarrange(plotlist = pp, ncol = length(pp), align = "h" )
dev.off()
## svg
## 2
ggpubr::ggarrange(plotlist = pp, ncol = length(pp), align = "h" )
p1 = plot.cluster.flow(mutdata,
cluster.col.name = cluster.col.name,
vaf.col.names = vaf.col.names,
sample.names =vaf.col.names,
colors = clone.colors,
y.title = "Variant Allele Frequency %"
) +
theme(axis.text.x = element_text(angle = 90) )
##
## Attaching package: 'reshape2'
## The following object is masked from 'package:tidyr':
##
## smiths
pdf( sprintf( "%s/%s.flow.pdf", output, output) , width=8, height=5 )
p1
dev.off()
## svg
## 2
p1
This script first infers the clonal orders in each sample by the pigeonhole principle. Then get the consensus trees across samples. Upon the built clonal trees, the scripts plot the clonal evolutionary history in each sample, the whole clonal structures for each sample, and the inferred clonal trees in the tumor.
source("data/auxiliary_functions.R")
# Inferring the main structures of mutation clusters
sample.groups = mapply(function(x) x[1], strsplit(vaf.col.names, "_") )
names(sample.groups) = vaf.col.names
sel = c(1:18)
y = plot.tree( output = output , mutdata = mutdata,
ccf.col.names = vaf.col.names[sel],
sample.groups = sample.groups[sel],
cancer.initiation.model = "monoclonal",
founding.cluster = 1 , ignore.clusters = 4 , cluster.col.name = "cluster",
subclonal.test.model = "non-parametric",
#scale.monoclonal.cell.frac = FALSE,
sum.p = 0.01, alpha = 0.05, weighted = FALSE,
plt.pairwise = F,
highlight.note.col.name = NULL,
highlight = "is.driver", highlight.CCF = TRUE
)
## infer.clonal.models
## Calculate VAF as CCF/2
## Sample 1: Coad_1 <-- Coad_1
## Sample 2: Coad_2 <-- Coad_2
## Sample 3: Coad_3 <-- Coad_3
## Sample 4: Coad_4 <-- Coad_4
## Sample 5: Coad_5 <-- Coad_5
## Sample 6: OveryLM_1 <-- OveryLM_1
## Sample 7: OveryLM_2 <-- OveryLM_2
## Sample 8: OveryLM_3 <-- OveryLM_3
## Sample 9: OveryLM_4 <-- OveryLM_4
## Sample 10: OveryLM_5 <-- OveryLM_5
## Sample 11: OveryRM_1 <-- OveryRM_1
## Sample 12: OveryRM_2 <-- OveryRM_2
## Sample 13: OveryRM_3 <-- OveryRM_3
## Sample 14: OveryRM_4 <-- OveryRM_4
## Sample 15: OveryRM_5 <-- OveryRM_5
## Sample 16: OveryRM_6 <-- OveryRM_6
## Sample 17: UterusM_1 <-- UterusM_1
## Sample 18: UterusM_3 <-- UterusM_3
## Using monoclonal model
## Note: all VAFs were divided by 100 to convert from percentage to proportion.
## Generating non-parametric boostrap samples...
## Coad_1 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 2,3,4,5,6,7,8
## User ignored clusters: 4
## Coad_1 : 2 clonal architecture model(s) found
##
## Coad_2 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 3,4,5,6,7,8
## User ignored clusters: 4
## Coad_2 : 6 clonal architecture model(s) found
##
## Coad_3 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 3,4,5,6,7,8
## User ignored clusters: 4
## Coad_3 : 6 clonal architecture model(s) found
##
## Coad_4 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 2,3,4,5,6,7,8
## User ignored clusters: 4
## Coad_4 : 2 clonal architecture model(s) found
##
## Coad_5 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 3,4,5,6,7,8,10
## User ignored clusters: 4
## Coad_5 : 2 clonal architecture model(s) found
##
## OveryLM_1 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,6,7,9,10
## User ignored clusters: 4
## OveryLM_1 : 12 clonal architecture model(s) found
##
## OveryLM_2 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,6,7,9,10
## User ignored clusters: 4
## OveryLM_2 : 10 clonal architecture model(s) found
##
## OveryLM_3 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,6,7,9,10
## User ignored clusters: 4
## OveryLM_3 : 13 clonal architecture model(s) found
##
## OveryLM_4 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,6,7,9,10
## User ignored clusters: 4
## OveryLM_4 : 12 clonal architecture model(s) found
##
## OveryLM_5 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,6,7,9,10
## User ignored clusters: 4
## OveryLM_5 : 13 clonal architecture model(s) found
##
## OveryRM_1 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,7,8,9,10
## User ignored clusters: 4
## OveryRM_1 : 13 clonal architecture model(s) found
##
## OveryRM_2 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,7,8,9,10
## User ignored clusters: 4
## OveryRM_2 : 13 clonal architecture model(s) found
##
## OveryRM_3 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,7,8,9,10
## User ignored clusters: 4
## OveryRM_3 : 12 clonal architecture model(s) found
##
## OveryRM_4 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,7,8,9,10
## User ignored clusters: 4
## OveryRM_4 : 10 clonal architecture model(s) found
##
## OveryRM_5 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,7,8,9,10
## User ignored clusters: 4
## OveryRM_5 : 12 clonal architecture model(s) found
##
## OveryRM_6 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 4,5,7,8,9,10
## User ignored clusters: 4
## OveryRM_6 : 12 clonal architecture model(s) found
##
## UterusM_1 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 2,3,5,6,7,8,9,10
## User ignored clusters: 4
## UterusM_1 : 1 clonal architecture model(s) found
##
## UterusM_3 : Enumerating clonal architectures...
## Determining if cluster VAF is significantly positive...
## Exluding clusters whose VAF < min.cluster.vaf=0.01
## Non-positive VAF clusters: 2,3,5,6,7,8,9,10
## User ignored clusters: 4
## UterusM_3 : 1 clonal architecture model(s) found
##
## Finding consensus models across samples...
## Found 6 consensus model(s)
## Generating consensus clonal evolution trees across samples...
## Found 6 consensus model(s)
## Scoring models...
## Pruning consensus clonal evolution trees....
## Seeding aware pruning is: off
## Number of unique pruned consensus trees: 6
## convert.consensus.tree.clone.to.branch
## plot.clonal.models
##
## Attaching package: 'igraph'
## The following objects are masked from 'package:dplyr':
##
## as_data_frame, groups, union
## The following objects are masked from 'package:purrr':
##
## compose, simplify
## The following object is masked from 'package:tidyr':
##
## crossing
## The following object is masked from 'package:tibble':
##
## as_data_frame
## The following objects are masked from 'package:stats':
##
## decompose, spectrum
## The following object is masked from 'package:base':
##
## union
## Plotting pruned consensus trees...
## Output plots are in: Met2
#PLot the clonal trees.
pdf(file = sprintf("%s/%s.trees.pdf", output, output), width = 6, height = 6)
plot.all.trees.clone.as.branch(y, branch.width = 0.5,
node.size = 2, node.label.size = 0.5,
tree.rotation =180,
angle = 20
)
dev.off()
## svg
## 2
plot.all.trees.clone.as.branch(y, branch.width = 0.5,
node.size = 2, node.label.size = 0.5,
tree.rotation =180,
angle = 20
)
plot.clonal.models(y,
# box plot parameters
box.plot = FALSE,
fancy.boxplot = FALSE,
# bell plot parameters
clone.shape = "bell",
bell.event = TRUE,
bell.event.label.color = "blue",
bell.event.label.angle = 60,
clone.time.step.scale = 1,
bell.curve.step = 2,
# node-based consensus tree parameters
merged.tree.plot = FALSE,
tree.node.label.split.character = NULL,
tree.node.shape = "circle",
tree.node.size = 30,
tree.node.text.size = 0.5,
merged.tree.node.size.scale = 1.25,
merged.tree.node.text.size.scale = 2.5,
merged.tree.cell.frac.ci = FALSE,
# branch-based consensus tree parameters
merged.tree.clone.as.branch = FALSE,
# cellular population parameters
cell.plot = TRUE,
num.cells = 100,
cell.border.size = 0.25,
cell.border.color = "black",
clone.grouping = "horizontal",
#meta-parameters
scale.monoclonal.cell.frac = TRUE,
show.score = FALSE,
cell.frac.ci = TRUE,
disable.cell.frac = FALSE,
# output figure parameters
out.dir = sprintf("%s/%s", output, "merge"),
out.format = "pdf",
overwrite.output = TRUE,
width = 9,
#height = ,
# vector of width scales for each panel from left to right
panel.widths = c(4,2))
## Plotting pruned consensus trees...
## Output plots are in: Met2/merge
save(mutdata, clone.colors, cluster.col.name , output, sampleNames, y, sample.groups, file = sprintf("%s/clonevol.data.rda", output) )
Merge multiple samples into a single sample. This can be used to merge multiregional sequencing (MRS) into one sample that could represent whole tumors.
Note: CCF should be in the range between 0-100
#merge coad
sel = 1:5
y.merge = merge.samples(y, samples = vaf.col.names[sel],
new.sample = "Coad", new.sample.group = "Coad",
ref.cols = str_c(vaf.col.names[sel], ".ref"),
var.cols = str_c(vaf.col.names[sel], ".var") )
## Note: all VAFs were divided by 100 to convert from percentage to proportion.
## Generating non-parametric boostrap samples...
## Estimating CCF of clones for merged sample...
## Merging clonal evolution trees across samples...
## WARN: Model 1 : Clone(s) 1 (existed in one of Coad_1+Coad_2+Coad_3+Coad_4+Coad_5 ) now lost in Coad
## ** merged sample Coad reannotated
## WARN: Model 2 : Clone(s) 1 (existed in one of Coad_1+Coad_2+Coad_3+Coad_4+Coad_5 ) now lost in Coad
## ** merged sample Coad reannotated
## WARN: Model 3 : Clone(s) 1 (existed in one of Coad_1+Coad_2+Coad_3+Coad_4+Coad_5 ) now lost in Coad
## ** merged sample Coad reannotated
## WARN: Model 4 : Clone(s) 1 (existed in one of Coad_1+Coad_2+Coad_3+Coad_4+Coad_5 ) now lost in Coad
## ** merged sample Coad reannotated
## WARN: Model 5 : Clone(s) 1 (existed in one of Coad_1+Coad_2+Coad_3+Coad_4+Coad_5 ) now lost in Coad
## ** merged sample Coad reannotated
## WARN: Model 6 : Clone(s) 1 (existed in one of Coad_1+Coad_2+Coad_3+Coad_4+Coad_5 ) now lost in Coad
## ** merged sample Coad reannotated
## Pruning consensus clonal evolution trees....
## Seeding aware pruning is: off
## Number of unique pruned consensus trees: 6
#merge overyLM
sel = 6:10
y.merge = merge.samples(y.merge, samples = vaf.col.names[sel],
new.sample = "OveryLM", new.sample.group = "OveryLM",
ref.cols = str_c(vaf.col.names[sel], ".ref"),
var.cols = str_c(vaf.col.names[sel], ".var") )
## Note: all VAFs were divided by 100 to convert from percentage to proportion.
## Generating non-parametric boostrap samples...
## Estimating CCF of clones for merged sample...
## Merging clonal evolution trees across samples...
## WARN: Model 1 : Clone(s) 2 (existed in one of OveryLM_1+OveryLM_2+OveryLM_3+OveryLM_4+OveryLM_5 ) now lost in OveryLM
## ** merged sample OveryLM reannotated
## WARN: Model 2 : Clone(s) 2 (existed in one of OveryLM_1+OveryLM_2+OveryLM_3+OveryLM_4+OveryLM_5 ) now lost in OveryLM
## ** merged sample OveryLM reannotated
## WARN: Model 3 : Clone(s) 2 (existed in one of OveryLM_1+OveryLM_2+OveryLM_3+OveryLM_4+OveryLM_5 ) now lost in OveryLM
## ** merged sample OveryLM reannotated
## WARN: Model 5 : Clone(s) 2 (existed in one of OveryLM_1+OveryLM_2+OveryLM_3+OveryLM_4+OveryLM_5 ) now lost in OveryLM
## ** merged sample OveryLM reannotated
## Pruning consensus clonal evolution trees....
## Seeding aware pruning is: off
## Number of unique pruned consensus trees: 3
#merge overyRM
sel = c(11:13, 14, 15, 16)
y.merge = merge.samples(y.merge, samples = vaf.col.names[sel],
new.sample = "OveryRM", new.sample.group = "OveryRM",
ref.cols = str_c(vaf.col.names[sel], ".ref"),
var.cols = str_c(vaf.col.names[sel], ".var") )
## Note: all VAFs were divided by 100 to convert from percentage to proportion.
## Generating non-parametric boostrap samples...
## Estimating CCF of clones for merged sample...
## Merging clonal evolution trees across samples...
## Pruning consensus clonal evolution trees....
## Seeding aware pruning is: off
## Number of unique pruned consensus trees: 2
#merge UterusM
sel = c(17, 18)
y.merge = merge.samples(y.merge, samples = vaf.col.names[sel],
new.sample = "UterusM", new.sample.group = "UterusM",
ref.cols = str_c(vaf.col.names[sel], ".ref"),
var.cols = str_c(vaf.col.names[sel], ".var") )
## Note: all VAFs were divided by 100 to convert from percentage to proportion.
## Generating non-parametric boostrap samples...
## Estimating CCF of clones for merged sample...
## Merging clonal evolution trees across samples...
## Pruning consensus clonal evolution trees....
## Seeding aware pruning is: off
## Number of unique pruned consensus trees: 2
plot.clonal.models(y.merge,
# box plot parameters
box.plot = FALSE,
fancy.boxplot = FALSE,
# bell plot parameters
clone.shape = "bell",
bell.event = TRUE,
bell.event.label.color = "blue",
bell.event.label.angle = 60,
clone.time.step.scale = 1,
bell.curve.step = 2,
# node-based consensus tree parameters
merged.tree.plot = TRUE,
tree.node.label.split.character = NULL,
tree.node.shape = "circle",
tree.node.size = 30,
tree.node.text.size = 0.5,
merged.tree.node.size.scale = 1.25,
merged.tree.node.text.size.scale = 2.5,
merged.tree.cell.frac.ci = FALSE,
# branch-based consensus tree parameters
merged.tree.clone.as.branch = FALSE,
# cellular population parameters
cell.plot = TRUE,
num.cells = 100,
cell.border.size = 0.25,
cell.border.color = "black",
clone.grouping = "horizontal",
#meta-parameters
scale.monoclonal.cell.frac = TRUE,
show.score = FALSE,
cell.frac.ci = TRUE,
disable.cell.frac = FALSE,
# output figure parameters
out.dir = sprintf("%s/%s", output, "merge1"),
out.format = "pdf",
overwrite.output = TRUE,
width = 9,
#height = ,
# vector of width scales for each panel from left to right
panel.widths = c(4,4,4))
## Plotting pruned consensus trees...
## Output plots are in: Met2/merge1
#merge MRS data into one.
save(y.merge , file = sprintf("%s/%s.y.merge.rda", output, output))